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Abstract 

A  main  challenge  in  transporting  cargo  for  United  States  Transportation  Com¬ 
mand  (USTRANSCOM)  is  in  mode  selection  or  integration.  Demand  for  cargo  is 
time  sensitive  and  must  be  fulfilled  by  an  established  due  date.  Since  these  clue  dates 
are  often  inflexible,  commercial  carriers  are  used  at  an  enormous  expense,  in  order 
to  fill  the  gap  in  organic  transportation  asset  capacity.  This  dissertation  develops 
a  new  methodology  for  transportation  capacity  assignment  to  routes  based  on  the 
Resource  Constrained  Shortest  Path  Problem  (RCSP).  Routes  can  be  single  or  mul¬ 
timodal  depending  on  the  characteristics  of  the  network,  delivery  timeline,  modal 
capacities,  and  costs.  The  difficulty  of  the  RCSP  requires  use  of  metaheuristics  to 
produce  solutions.  An  Ant  Colony  System  to  solve  the  RCSP  is  developed  in  this 
dissertation.  Finally,  a  method  for  generating  near  Pareto  optimal  solutions  with 
respect  to  the  objectives  of  cost  and  time  is  developed. 
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Rough-Cut  Capacity  Planning  in  Multimodal  Freight 


Transportation  Networks 


I.  Introduction 

1.1  Background,  &  Motivation 

Transportation  networks  are  prevalent  in  nearly  every  aspect  of  our  society 
and  as  of  the  most  recent  economic  census,  transportation  revenues  in  the  U.S.  rep¬ 
resent  a  nearly  $640  billion  segment  of  the  economy  [1],  Domestic  and  international 
commerce  is  dependent  upon  the  transportation  of  goods,  personnel,  and  commu¬ 
nications  between  geographically  diverse  locations.  Commercial  airlines,  freight  ser¬ 
vices,  and  internet-based  auction  sites  all  depend  on  access  to  efficient  and  effective 
transportation  networks  to  maintain  a  competitive  edge  and  continue  profitable  op¬ 
erations. 

Profitability  of  civilian  industry  depends  heavily  on  the  agility  and  cost  effec¬ 
tiveness  of  the  transportation  and  supply  chains  used  to  transport  raw  materials 
for  production  and  to  ship  finished  goods  to  market.  This  need  has  led  to  the 
establishment  and  growth  of  Third  Party  Logistics  providers  (3PL).  These  3PL  or¬ 
ganizations  and  DoD  supply  chain  counterparts  are  under  continuous  pressure  to 
provide  increasingly  efficient  transportation  options  to  their  customers.  This  can  be 
attributed  to  globalization,  shrinking  budgets,  and  expanding  mission  requirements. 
Meeting  these  demands  requires  distribution  providers  to  integrate  historically  dis¬ 
parate  transportation  modes  creating  true  multimodal  options. 

Military  applications  of  transportation  are  equally  dependent  on  the  perfor¬ 
mance  of  the  networks  they  utilize.  The  problem  of  transporting  commodities,  and 
people  is  so  important  to  the  DoD,  that  the  DoD  has  attempted  to  centralize  trans- 
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portation  control  across  the  Army,  Navy,  and  Air  Force  as  early  as  1979  with  the 
formation  of  the  Joint  Deployment  Agency  [90].  In  1987  the  United  States  Trans¬ 
portation  Command  (USTRANSCOM)  was  founded.  Today  the  USTRANSCOM 
is  responsible  for  all  assets  and  commodities  via  the  Surface  Deployment  and  Dis¬ 
tribution  Command  (SDDC),  Military  Sealift  Command  (MSC),  and  Air  Mobility 
Command  (AMC)  until  it  reaches  an  Area  of  Responsibility  (AOR)  at  which  time 
responsibility  for  goods  in  transit  is  assumed  by  a  Geographic  Combatant  Comman¬ 
der  (GCC).  In  order  to  support  these  operations,  USTRANSCOM  must  decide  how 
to  allocate  the  air,  land  (truck/rail),  and  sea  assets  as  well  as  personnel  resources  in 
such  a  way  as  to  ensure  cargo  is  delivered  on  time  to  the  locations  it  services  [81].  The 
2011  USTRANSCOM  Strategic  Plan  discusses  the  need  for  JDDE  synchronization 
across  all  modes  of  transportation  utilized  by  the  JDDE  partners  [107]. 

Although  multimodal  transportation  networks  are  used  every  day  in  both  mil¬ 
itary  and  civilian  applications,  the  process  of  planning  and  coordinating  transporta¬ 
tion  modes  has  received  surprisingly  little  attention  in  the  literature.  Many  recent 
publications  in  the  held  of  Supply  Chain  Management  and  Transportation  claim  to 
address  multimodal  transportation  planning  but,  in  fact  only  address  the  coordina¬ 
tion  of  at  most  two  transportation  modes  [52],  Addressing  multimodal  transporta¬ 
tion  as  an  integrated  network  system  is  a  relatively  new  concept.  In  fact  the  first 
instance  of  this  systems  approach  found  in  the  literature  is  from  1999  [88].  This 
research  addresses  capacity  planning  in  multimodal  freight  networks. 

1.2  Problem,  Description,  Themes,  &  Gaps 

Planning  and  coordinating  multimodal  transportation  networks  requires  a  method 
of  estimating  the  ability  of  the  network  to  meet  demands.  This  estimate  must  ex¬ 
plicitly  consider  many  factors  affecting  the  overall  capacity  that  the  network  can 
support.  Network  capacity  over  a  specified  time  period  is  sometimes  referred  to  as 
throughput.  This  definition  of  throughput  will  be  adopted  in  the  remainder  of  this 
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research.  Freight  network  throughput  is  affected  by  many  factors.  Those  considered 
in  this  research  are:  planning  horizon  length,  fleet  size  for  each  transportation  mode, 
payload  (in  tons)  of  each  asset  by  type,  asset  utilization  rate,  block  speed,  and  asset 
productivity. 

Not  only  is  the  problem  of  estimating  capacity  in  transportation  networks 
difficult,  but  given  the  increased  demands  for  network  capacity  dne  to  increases  in 
cargo  volume  demand,  the  problem  has  become  a  major  concern  of  public  policy  [88]. 

Estimating  the  capacity  of  a  multimodal  network  is  plagued  by  the  inadequa¬ 
cies  in  current  capacity  modeling  methods  [96].  Previous  work  in  estimating  the 
capacity  of  multimodal  freight  networks  is  inaccurate  due  to  measuring  individual 
system  components  or  the  capacity  of  partially  aggregated  components.  These  par¬ 
tial  aggregations  usually  involve  aggregation  for  estimating  the  capacity  of  single 
transportation  modes  or  aggregation  in  estimating  the  capacity  of  single  “corridors” 
through  the  multimodal  system  [49,88,105,115].  The  inadequacy  stems  from  the 
inability  of  current  models  to  capure  interactions  among  transportation  modes  due 
to  mid-path  mode  transfers  and  the  interaction  of  other  system  components  [96]. 

The  most  recent  improvements  in  the  area  of  multimodal  freight  transporta¬ 
tion  system  capacity  seek  to  provide  a  systems  view  of  multimodal  transportation 
networks  but  addresses  the  problem  at  a  tactical  level  of  detail  and  only  provides 
a  model  formulation.  No  attempt  is  made  to  develop  a  solution  approach  for  the 
bilevel  programming  problem  resulting  from  the  model  formulation  [96]. 

1.2.1  Themes.  Three  themes  provide  context  for  the  existing  literature 
and  contributions:  Integration  of  multiple  fleets  of  transportation  across  multiple 
modes  of  transportation  and  mode/carrier  selection.  Modal  integration  addresses 
the  problem  of  coordinating  the  resources  of  the  different  transportation  modes  in 
order  to  gain  efficiencies  through  synergy.  The  literature  review  will  demonstrate 
that  modal  integration  is  typically  limited  to  two  modes  of  transportation.  Mode 
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selection  seeks  the  single  “best”  mode  of  transportation  with  respect  to  cost,  security, 
delivery  time,  or  some  combination  of  these  and  other  objectives. 

Another  emergent  theme  that  is  not  addressed  in  the  literature  in  the  context  of 
multimodal  transportation  is  the  identification  of  efficient  or  near-efficient  solutions 
to  capacity  problems.  An  efficient  or  Pareto  optimal  solution  is  a  multi-objective 
concept  of  optimality.  An  efficient  solution  is  one  which  cannot  be  improved  with 
respect  to  one  objective  without  sacrificing  quality  in  at  least  one  other  objective  [40]. 
Efficiency  will  be  discussed  further  in  Chapter  5. 

1.2.2  Gaps.  In  the  held  of  production  and  inventory  management,  there  is 
a  keen  sense  of  the  fidelity  required  to  execute  planning  during  the  long-term,  mid¬ 
term,  and  short-term  planning  horizons.  This  awareness  helps  to  control  extraneous 
data  collections  which  are  often  both  costly  and  time-consuming.  Data  for  high 
fidelity  planning  can  be  prohibitively  difficult  to  obtain.  Additionally,  high  fidelity 
planning  (scheduling  the  production  on  a  single  machine  over  the  course  of  an  hour 
for  instance),  is  not  required  in  order  to  plan  production  operations  over  the  course 
of  the  next  year.  The  purpose  of  the  long  and  mid-term  planning  is  to  estimate 
the  feasibility  of  a  desired  level  of  production  relative  to  the  manufacturing  system 
capacity.  Long-term  (aggregate)  plans  are  refined  via  mid-term  (rough-cut  capacity) 
planning  as  the  objective  planning  period  draws  closer.  When  the  objective  planning 
period  is  within  a  month,  a  master  production  schedule  is  created  from  the  rough- 
cut  estimates  obtained  through  mid-term  and  long-term  planning.  In  this  way,  the 
production  capacity  of  the  system  is  not  overwhelmed  in  any  period. 

No  current  source  in  the  literature  provides  a  time  scalable  method  of  rough- 
cut  multimodal  freight  transportation  network  capacity  estimation.  Scalability  refers 
to  the  ability  of  the  method  to  adapt  capacity  estimates  to  an  arbitrarily  defined 
planning  horizon.  Furthermore,  the  limited  number  of  capacity  estimation  methods 
are  based  on  tactical-level  considerations  and  are  not  well  suited  to  estimating  overall 
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system  capacity  without  tactical-level  data  which  may  not  be  available  [96].  This 
observation  led  to  the  identification  of  several  gaps  in  the  literature  addressed  by 
this  research: 

•  (Gl)  No  scalable  multimodal  freight  transportation  network  capacity  estima¬ 
tion  method  is  found  in  the  currently  available  literature. 

•  (G2)  There  is  a  lack  of  an  operational-level  multimodal  freight  transporta¬ 
tion  capacity  estimation  method  requiring  only  operational-level  data.  Such  a 
method  would  be  an  extension  of  “capacity  planning  based  on  overall  factors” 
(CPOF)  used  in  production  and  inventory  management. 

•  (G3)  Currently,  there  is  no  rough-cut  capacity  planning  technique  for  mul¬ 
timodal  transportation.  Developing  such  a  capability  by  extending  existing 
concepts  would  represent  a  major  contribution  to  the  field  of  transportation 
research. 

•  (G4)  Development  of  a  new  modeling  approach  that  considers  a  multimodal 
transportation  system  holistically  would  represent  a  contribution  to  the  field 
of  transportation  research. 

•  (G5)  Development  of  a  solution  procedure  to  identify  efficient  solutions  under 
various  planning  assumption  scenarios  would  represent  a  contribution  to  the 
fields  of  operations  research  and  transportation  research. 

1.3  Problem  Statement 

Identification  of  these  gaps  leads  immediately  to  the  following  problem  state¬ 
ment:  “Given  a  fleet  of  vehicles  from  across  multiple  modes  of  transportation,  de¬ 
mand  quantities/locations,  supply  quantities/locations,  and  planning  horizon  length, 
determine  the  multimodal  fleet  mix  and  route  providing  the  required  capacity  at 
minimum  cost  based  only  upon  these  operational-level  data  inputs” 
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The  capacity  of  a  network  relates  to  its  throughput  over  a  specified  period  of 
time.  Since  network  flow  may  be  continuous,  the  amount  of  throughput  must  be 
defined  for  a  specified  time  period  in  order  to  provide  a  basis  for  measurement  and 
comparison.  Capacity  estimates  in  freight  transportation  networks  are  conveniently 
expressed  in  units  of  million  ton  miles  per  day  (MTM/D).  This  standard  measure  of 
vehicle  fleet  performance  is  the  standard  for  air  mobility  but  may  be  calculated  for 
any  fleet  of  vehicles  [91]. 

1.3.1  Research  Questions.  Effectively  addressing  this  problem  statement 
leads  to  several  questions  that  must  be  answered.  The  specific  questions  of  interest 
in  this  research  are: 

•  (Q1.G1/G2/G3)  Can  concepts  of  Rough-Cut  Capacity  Planning  be  extended 
to  estimate  the  capacity  of  a  multimodal  freight  transportation  network  as  a 
system  given  existing  routes? 

•  (Q2.G2/G3)  Can  modeling  and  solution  techniques  be  developed  to  allocate 
the  available  capacity  of  a  multimodal  freight  transportation  network  in  a  way 
that  minimizes  costs  and  meets  demand  timelines? 

•  (Q3.G2/G3)  What  multimodal  routes  should  be  established  in  a  network,  in 
order  to  improve  cost  efficiency  and  mission  effectiveness? 

•  (Q4.G5)  How  do  changes  in  numbers  of  assets,  available  routes,  planning  hori¬ 
zon  length,  and  input  costs  (fuel  for  instance)  affect  mode  choices  and  overall 
system  costs? 

•  (Q5.G4)  What  areas  of  the  multimodal  transportation  system  would  benefit 
most  from  resource  investments? 

•  (Q6.G1/G2)  What  is  the  appropriate  fieet  mix  required  to  fulfill  demand  in  a 
set  planning  horizon  at  a  minimum  cost  to  the  overall  system? 
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•  (Q7.G5)  How  robust  is  the  multimodal  transportation  system  to  losses  in  es¬ 
tablished  routes  or  assets?  That  is,  what  is  the  nature  and  length  of  system 
capability  degradation  due  to  these  losses? 

1-4  Research  Objectives 

The  hypothesis  of  this  research  is  that  multimodal  rough-cut  system  capacity 
planning  can  be  executed  using  a  network-based  approach  and  capacity  planning 
techniques  normally  applied  to  production  and  inventory  management.  Extending 
these  methods  requires  the  development  of  new  theory.  The  following  research  objec¬ 
tives  provide  original  theoretical  contributions  and  move  the  state-of-the-art  forward 
with  respect  to  closing  the  gaps  identified  above. 

Addressing  these  questions  requires  the  development  of  a  new  methodology  and 
new  modeling  and  analysis  approaches  to  augment  what  is  currently  available  in  the 
literature.  These  questions  lead  to  an  explicit  statement  of  the  research  objectives. 

I.4.I  Research  Objectives  &  Contributions. 

•  (01. Q1  / Q2/Q3)  Develop  a  methodology  for  operational-level  multimodal  freight 
transportation  capacity  estimation  that  is  scalable  to  an  arbitrary  planing  hori¬ 
zon  and  an  arbitrary  multimodal  fleet  size.  The  modeling  method  should  cap¬ 
ture  interaction  between  the  various  modes  of  transportation. 

•  (02.Q4/Q5/Q6)  Develop  a  method  for  determining  the  least  cost  fleet  mix 
and  routing  structure  for  a  specified  range  of  assets  and  a  specified  range  on 
planning  horizon. 

•  (03. Q7)  Develop  a  method  for  generating  alternative  near  efficient  solutions 
and  solutions  with  passive  resilience  to  aid  decision  makers  in  analyzing  trade¬ 
offs  in  time,  cost,  and  fleet  assets. 
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1.5  Dissertation  Overview 


Chapter  1  has  motivated  this  dissertation  effort  by  briefly  highlighting  the  na¬ 
ture  of  the  problem,  discussing  the  relevance  of  the  research  within  the  context  of 
both  military  and  commercial  applications,  and  has  presented  the  intended  contri¬ 
butions. 

Chapter  2  will  establish  the  uniqueness  of  the  research  in  this  dissertation  by 
providing  a  thorough  treatment  of  the  relevant  literature.  This  chapter  will  compare 
and  contrast  this  dissertation  work  with  past  work.  This  chapter  shows  that  the 
intended  contributions  are  indeed  original  and  provide  a  solid  scope  and  context  for 
the  Methodology. 

In  chapter  3  we  show  how  Multimodal  Rough-Cut  Capacity  Planning  is  mod¬ 
eled  using  the  Resource  Constrained  Shortest  Path  Problem.  We  demonstrate  how 
this  approach  supports  either  mode  selection  or  integration  depending  upon  trans¬ 
portation  mode  costs,  fleet  size,  and  planning  horizon. 

In  chapter  4  we  develop  an  Ant  Colony  System  metaheuristic  to  quickly  solve 
large  instances  of  the  RCSP.  In  chapter  5  we  develop  a  heuristic  based  on  the  ACS 
of  chapter  4.  This  new  metaheuristic  generates  near-efficient  solutions  and  helps 
identify  passive  resilient  solutions  to  the  multimodal  RCCP  problem  to  aid  decision 
makers  in  analyzing  tradeoffs  between  time,  cost  and  asset  utilization. 

Chapter  6  summarizes  the  contributions,  results,  and  conclusions  of  the  disser¬ 
tation  and  presents  final  analyses.  Recommendations  for  future  work  complete  this 
chapter. 


II.  Literature  Review 


2.1  Literature  Review  Strategy 

This  chapter  provides  a  broad  overview  of  the  relevant  literature.  Chapters 
3,  4,  and  5  contain  context  specific  literature  review  and  thus  the  purpose  of  this 
chapter  is  to  demonstrate  the  uniqueness  of  the  dissertation  research  described  in 
Chapter  1. 

2.2  Overview  of  the  Literature 

Figures  2.1  and  2.2  provide  a  quick  reference  summary  of  the  important  sources 
reviewed  in  this  dissertation.  Some  sources,  such  as  textbooks,  were  omitted  since 
their  purpose  is  to  support  the  analysis  rather  than  contribute  original  work  to  any 
held  of  research. 
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Multi- 
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Time/- 
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Math 

Programming 

X 

X 

X 
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X 

Scheduling 

Heuristics 

X 

Tabic  2.1:  Features  of  this  dissertation  versus  literature’s  research 


Capacity  planning  in  transportation  is  usually  limited  to  one  or  two  modes 
of  transportation.  Sources  considering  more  than  two  modes  of  transportation  are 
considering  the  problem  of  mode  selection  [83,96].  True  multimodal  integration  is  a 
topic  which  is  considered  an  area  of  opportunity  in  recent  surveys  of  the  transporta¬ 
tion  field  [79,119].  The  focus  of  this  dissertation  is  on  the  development  of  methods 
for  capacity  analysis  which  are  scalable  to  an  arbitrary  number  of  transportation 
modes.  Table  2.1  provides  a  concise  comparison  of  this  dissertation  with  other  areas 
of  reviewed  literature. 

Sources  in  the  literature  focus  on  either  mode  selection  or  mode  integration. 
This  dissertation  is  unique  in  several  ways  including  the  capability  of  analyzing  a 
multimodal  network  considering  both  mode  selection  and  mode  integration.  A  single 
mode  solution  is  returned  if  it  is  the  most  cost  effective  or  is  “best”  with  respect  to 
the  multiple  defined  objectives.  A  multimodal  solution  is  returned  if  it  is  the  most 
cost  effective,  is  the  best  with  respect  to  the  multiple  objectives,  or  if  multiple  modes 
of  transportation  are  needed  in  order  to  make  the  problem  feasible  by  providing  the 
required  capacity  to  meet  delivery  demands  in  the  specified  planning  horizon.  Table 
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Chapter 

Unique 

Transportation 
Planning  Capability 

Ops  Research 
Methodology 

3 

Multimodal  RCCP 

RCSP  Problem 

4 

Accessibility  of  Large 
Multimodal  RCCP 
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5 

Generate  Alternative 
Near- Efficient  and 
Passive  Resilient 
Solutions 

New  Heuristic  based  on 
ACS 

Tabic  2.2:  Overview  of  this  dissertation 


2.2  presents  a  summary  of  the  content  of  the  dissertation  development  in  chapters 
3,  4,  and  5. 

Capacity  planning  can  be  viewed  from  either  the  aggregate  or  the  rough-cut 
perspective.  Aggregate  planning  focuses  on  the  longer  term  strategy  of  production 
and  attempts  to  ensure  supply  matches  demand  by  prescribing  monthly  or  weekly 
production  output  levels.  By  contrast,  rough-cut  planning  seeks  to  ensure  that  the 
production  level  prescribed  in  the  aggregate  planning  process  is  feasible  by  estimating 
the  potential  throughput  or  capacity  of  the  production  system. 


2.3  Multimodal  Freight  Transportation 

Literature  sources  reviewed  in  this  research  use  the  terms  intermodal  and  mul¬ 
timodal  interchangeably.  In  order  to  standardize  their  use  in  this  dissertation  the 
terms  shall  be  defined  as  follows: 


•  Multimodal  transportation  refers  to  the  delivery  of  freight  from  origin  to  des¬ 
tination  involving  at  least  two  different  modes  of  transportation. 
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•  Intermodal  transportation  is  refers  to  the  transfer  operations  to  move  freight 
between  two  modes  of  transportation  at  container  terminals  or  other  types  of 
intermodal  ports  [23]. 

2.3.1  Intermodal  Transportation.  Intermodal  transportation  is  defined  as: 
“the  transportation  of  a  person  or  load  from  its  origin  to  its  destination  by  a  se¬ 
quence  of  at  least  two  transportation  modes,  the  transfer  from  one  mode  to  the  next 
being  performed  at  an  intermodal  terminal”  [29].  The  process  of  moving  the  goods 
between  modes  is  often  referred  to  as  transshipment.  It  is  not  a  goal  of  the  research 
which  follows  to  treat  transshipment  with  any  level  of  detail.  Rather,  this  section  is 
designed  to  acknowledge  that  this  topic  is  indeed  an  area  of  research  in  itself.  For  fur¬ 
ther  review  of  the  topic  of  transshipment  the  reader  is  referred  to  [9,76,85,108,112], 
These  efforts  are  targeted  at  the  tactical  level  and  focus  on  port  operations  including 
how  to  utilize  different  technologies  to  efficiently  accomplish  port  tasks.  For  the  pur¬ 
poses  of  this  research,  tactical  efficiencies  matter  only  as  they  relate  to  throughput 
capacities  of  the  various  nodes  (bases)  in  the  network. 

A  1996  article  by  Haghani  and  Sei-Chang  models  logistics  for  disaster  relief  as 
a  multi-commodity,  multi-modal  network  flow  problem  [53] .  The  modeling  approach 
utilizes  the  concept  of  a  time-space  network  which  has  also  been  used  by  Nielsen  et 
al.  and  in  many  other  applications  by  Barnhart  and  Armacost  in  a  similar  manner  [6, 
7,95].  The  authors  formulate  the  network  as  a  MILP  and  then  propose  two  different 
heuristics  for  obtaining  solutions  to  the  formulation.  Route  length  restrictions  are 
not  considered  in  this  formulation.  By  using  the  time-space  network  formulation, 
the  model  is  capable  of  representing  a  change  in  transportation  mode  by  an  arc 
traversal.  Similarly,  advances  in  time  are  represented  as  arc  traversals  as  well.  The 
distinction  between  changes  in  mode,  changes  in  time,  and  changes  in  location  are 
managed  by  indexing  the  sets  of  time-windows,  locations,  and  modes. 
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2-4  Multi-Layered  Networks 

Kennedy  introduced  the  concept  of  multi-layered  network  synthesis  considering 
both  cost  and  robustness  as  competing  objectives.  He  applied  a  multicriteria  decision 
analysis  technique  using  epsilon  constraints  to  scalarize  the  objective  function.  He 
then  extended  this  new  formulation  to  include  flow  and  connectivity  constraints  [73]. 
Multimodal  freight  transportation  networks  can  be  modeled  using  layered  networks 
and  interactions  (i.e.  transfers  between  modes)  can  be  modeled  using  inter-layer 
arcs  [73,74], 

2.5  Hierarchical  Production  Planning 

Hierarchical  Production  Planning  techniques  have  been  developed  for  use  in 
planning  production  operations  at  varying  levels  of  fidelity.  Long-range  planning 
is  generally  executed  using  lower  levels  of  fidelity  since  there  may  be  changes  in 
important  factors  affecting  production  between  the  planning  and  execution  phases. 
Mid-range  planning  is  somewhat  higher  fidelity  but  still  does  not  produce  detailed 
tactical  level  schedules  due  to  potential  changes  in  input  factors.  Master  Production 
Scheduling  produces  detailed  schedules  of  how  the  raw  inputs  will  be  converted  into 
finished  products  using  the  manufacturing  resources  of  the  factory,  like  machines  and 
manpower  [19, 30, 58, 60, 61,69, 93, 98]. 

2.5.1  Aggregate  Production  Planning.  Aggregate  planning  encompasses 
both  aggregate  output  planning  and  aggregate  capacity  planning.  Aggregate  output 
planning  attempts  to  match  production  and  demand.  Aggregate  capacity  planning 
assesses  the  feasibility  of  the  aggregate  output  plan  by  estimating  how  much  produc¬ 
tion  capacity  is  required  to  meet  the  production  levels  prescribed  by  the  aggregate 
output  plan  [69] .  In  aggregate  planning  the  individual  products  are  normally  aggre¬ 
gated  into  product  families  based  on  some  measure  of  product  similarity.  By  reducing 
the  number  of  products  the  aggregate  planning  process  can  be  completed  in  a  timely 
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manner  [28].  The  estimates  obtained  using  aggregate  planning  are  rough,  but  are 
usually  sufficient  to  estimate  the  long-term  requirements  in  the  six  to  eighteen  month 
rolling  planning  horizon. 

2.5.2  Rough-Cut  Capacity  Planning.  The  focus  of  this  dissertation  is 
rough-cut  capacity  planning  in  multimodal  transportation.  Rough-cut  capacity  plan¬ 
ning  is  similar  to  aggregate  capacity  planning  but  differs  in  fidelity  and  function. 
In  rough-cut  planning  the  goal  is  to  estimate  the  ability  of  the  system  to  meet 
the  required  production  demand  over  a  shorter  time  period  and  is  used  to  assess 
the  feasibility  of  a  master  production  schedule.  RCCP  addresses  the  question:  “Is 
there  sufficient  capacity  in  the  production  capabilities  to  support  the  desired  level 
of  production?”  There  are  three  prominent  methods  of  performing  RCCP:  Capacity 
Planning  using  Overall  Factors  (CPOF),  Capacity  Bill  Procedures,  Resource  Pro¬ 
files  [66,109], 

In  the  context  of  this  dissertation,  rough-cut  capacity  planning  will  be  de¬ 
fined  as  planning  to  establish  the  feasibility  of  a  desired  level  of  multimodal  system 
freight  throughput  in  millions  of  tons.  The  scheduling  of  the  multimodal  fleet  is 
not  addressed  but  many  heuristics  for  scheduling  fleets  exist  and  the  reader  is  en¬ 
couraged  to  reference  Pinedo  for  a  thorough  treatment  of  scheduling  algorithms  and 
heuristics  [97]. 

The  function  of  the  capacity  planning  and  the  data  requirements  are  the  pri¬ 
mary  reasons  for  classifying  this  dissertation  as  rough-cut  capacity  planning  rather 
than  aggregate  capacity  planning.  As  Chapter  3  will  show,  only  operational  level 
data  such  as  aircraft  capacity,  and  speed  are  needed  in  order  to  complete  the  capac¬ 
ity  planning  approach  developed  in  this  research.  Also,  the  heuristic  developed  is 
meant  to  provide  a  user  with  the  ability  to  quickly  produce  solutions  under  differ¬ 
ent  scenarios  so  the  function  is  more  in  line  with  near-term  rather  than  long-term 
planning.  However,  the  math  programming  model  provided  in  chapter  3  may  also 
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be  solved  by  deterministic  methods  if  time  allows  so  its  scope  is  not  restricted  only 
to  near  term  estimates  as  optimal  solutions  may  be  determined  as  well. 

2.5.3  Master  Production  Scheduling.  The  purpose  of  master  production 
scheduling  is  to  plan  the  short-term  (i.e.,  weekly  or  monthly)  production  activities 
of  the  factory  in  order  to  ensure  demand  is  met.  In  MPS  the  product  families  are  dis¬ 
aggregated  into  individual  products  and  are  treated  individually.  A  key  requirement 
for  the  successful  linking  of  aggregate  planning,  rough-cut  planning,  and  master  pro¬ 
duction  scheduling  is  in  aggregating  and  disaggregating  products  and  product  fami¬ 
lies  in  an  appropriate  way.  For  more  detail  on  the  issues  associated  with  aggregation 
and  disaggregation  the  reader  is  encouraged  to  reference  Hax  and  Candea  [58]. 

2.6  Capacity  Planning  in  Multimodal  Freight  Transportation 

Park  developed  a  hrst-of-its-kind  model  to  address  multimodal  freight  trans¬ 
portation  from  a  comprehensive  standpoint  in  order  to  leverage  the  competitive 
advantages  of  each  mode  [96].  The  model  Park  develops  is  an  integrated  multimodal 
freight  model  but  the  model  uses  an  “average  load  factor”  for  each  mode  and  there¬ 
fore  is  not  additive  with  respect  to  capacity  over  the  entire  fleet  of  vehicles  in  each 
mode.  This  is  in  contrast  to  the  MTM/D  approach  proposed  in  this  research  [96]. 

Sun  et.al.  use  an  existing  LP  model  called  MAXCAP  which  seeks  to  maximize 
the  system  capacity  for  multiple  origin  destination  pairs  subject  to  resource  and 
capacity  constraints.  Additionally,  they  analyze  the  “flexibility”  and  performance 
of  the  system  under  degradation.  However,  they  consider  only  a  single  mode  of 
transportation,  movement  by  rail  [104], 

Andersen  et.  al.  address  the  design  and  optimization  of  a  multimodal  sched¬ 
uled  service  network.  The  model  provides  “tactical”  level  planning.  The  main  goal 
of  the  model  developed  is  to  determine  the  appropriate  departure  times  to  maximize 
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throughput.  The  primary  issue  in  this  report  is  the  coordination  and  synchronization 
of  schedules  among  multiple  fleets  of  vehicles  [4], 

Kasikitwiwat  and  Chen  discuss  different  definitions  for  network  capacity.  The 
authors  categorize  capacity  models  into  two  main  types:  economic  and  physical  net¬ 
work  capacity  [71].  Three  methods  of  examining  physical  network  capacity  (reserve, 
ultimate,  and  practical  capacities)  are  presented  and  discussed. 

Russ  et.al.  formulate  a  network  design  problem  for  multimodal  transportation 
infrastructure  investment  planning  in  Indonesia,  as  a  bi-level  programming  problem. 
The  lower-level  program  models  traffic  flow  equilibrium  and  is  solved  using  a  iter¬ 
ative  descent  direction  algorithm.  The  upper  level  program  manages  the  capacity 
allocation  on  the  arcs  and  is  solved  using  a  genetic  algorithm  [102], 

2.1  Shortest  Path  and  Resource  Constrained  Shortest  Path  Problems 

The  RCSP  is  a  generalization  of  the  well  known  shortest  path  problem  in 
which  the  elements  of  the  shortest  path  are  also  subject  to  knapsack  constraints 
based  upon  resource  limitations.  The  RCSP  was  first  proven  to  be  NP-Complete  in 
1980  by  Handler  and  Zang  and  many  extensions  and  applications  have  since  been 
developed  [38,48,54],  Even  the  least  general  version  of  the  RCSP,  the  Elementary 
RCSP  (ERCSP)  in  which  no  arc  is  repeated  has  been  shown  to  be  strongly  NP-Hard 
by  Beasley  and  Christohdes  in  1989  [18].  This  ERCSP  Ends  applications  in  vehicle 
routing  where  it  is  often  the  route  pricing  subproblcm  [26,51,101]. 
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III.  Rough- Cut  Capacity  Planning  using  Resource- Constrained 

Shortest  Paths 

3. 1  Introduction 

Consider  the  problem  of  planning  for  capacity  utilization  in  multimodal  freight 
distribution  networks.  This  complex  problem  has  received  relatively  little  attention 
in  the  literature  from  a  truly  integrated  multimodal  perspective  [79].  Capacity  alloca¬ 
tion  in  multimodal  freight  transportation  must  address  many  factors  in  determining 
to  which  routes  capacity  is  to  be  allocated. 

The  literature  sources  reviewed  demonstrate  a  gap  in  the  area  of  transportation 
mode  integration.  Much  of  the  previous  research  addressing  multiple  transportation 
modes  focuses  on  mode  selection  rather  than  integration.  Mode  integration  is  defined 
here  as  the  seamless  use  of  two  or  more  modes  of  transportation  to  provide  required 
system  capacity.  Mode  integration  differs  from  mode  selection  in  that  mode  selection 
is  concerned  with  what  percentage  of  demand  is  to  be  allocated  to  each  mode  of 
transportation,  usually  based  upon  cost.  These  decision  problems  are  typical  in 
third  party  logistics  (3PL)  applications.  A  weakness  of  this  approach  is  that  mode 
selection  views  transportation  modes  in  a  strictly  parallel  configuration.  This  view  of 
the  transportation  modes  as  separate  systems  ignores  the  possibility  of  cooperation 
among  multiple  modes  to  save  time  and  money. 

A  recent  survey  of  the  transportation  mode  choice  and  carrier  selection  lit¬ 
erature  reviewed  over  ninety  separate  references  spanning  a  period  of  over  forty 
years  [83].  These  efforts  seek  the  mode  that  provides  the  least  expensive  shipping 
option  to  support  delivery  of  the  required  freight  capacity  within  the  specified  time 
period.  Research  efforts  following  this  methodology  ignore  the  potential  for  even 
greater  cost  and  time  savings  that  might  be  gained  through  truly  integrated  multi¬ 
modal  shipping  options. 
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Many  articles  have  been  devoted  to  estimating  the  capacity  of  single  mode 
transportation  networks  and,  in  particular,  rail  network  capacity  has  been  exten¬ 
sively  studied.  Many  of  these  efforts  use  or  extend  the  MAXCAP  model  which  uses 
“traffic  lanes”  that  are  aggregated  network  nodes  geographically  close  with  an  estab¬ 
lished  recurring  commodity  flow  [87,88,104],  Some  “multimodal”  capacity  planning 
sources  surveyed  were,  in  fact,  only  bi- modal  [16].  This  is  an  obvious  inadequacy  in 
the  literature  since  many  logistics  providers  use  a  combination  of  air,  truck,  rail,  and 
sea  assets.  Methods  of  simultaneous  mode  analysis  are  badly  needed. 

Previous  network  capacity  modeling  efforts  are  inadequate  to  capture  the  true 
transportation  network  capacity  since  they  typically  utilize  network  aggregation 
techniques  which  reduce  the  fidelity  and  accuracy  of  the  resulting  capacity  esti¬ 
mate.  Corridors  and  modal  networks  are  two  popular  methods  of  aggregation  used 
in  estimating  single  and  multimodal  network  capacities  [49,88,105].  A  recent  model¬ 
ing  effort  focusing  on  multimodal  transportation  capacity  estimation  uses  a  bi-level 
programming  formulation  to  allocate  capacity.  The  formulation  is  so  complex  no 
attempt  is  made  to  provide  a  solution  method  [96].  Aggregation  is  done  a  priori 
resulting  in  a  network  abstraction  which  has  no  ability  to  yield  insights  at  higher 
fidelity.  The  final  result  is  an  inadequate  analysis  resulting  from  a  loss  of  valuable 
system  information  due  to  aggregation. 

Macharis  and  Bontekoning  have  established  opportunities  in  operations  re¬ 
search  for  contributions  to  multimodal  transportation  [79].  The  gaps  in  the  litera¬ 
ture  are  categorized  into  drayage,  terminal,  network,  and  intermodal  categories  based 
upon  the  level  of  problem  fidelity.  This  work  considers  the  “network”  and  “inter¬ 
modal”  areas.  Drayage  is  mainly  used  in  reference  to  the  transport  of  containerized 
cargo  over  short  distances.  Terminal  operations  are  the  port  operations  concerned 
with  transloading  cargo  between  modes  or  vessels  of  conveyance  Network  issues  con¬ 
cern  infrastructure  and  route  planning  and  pricing.  Intermodal  or  multimodal  issues 
surround  the  route  and  service  choices  for  existing  networks.  In  particular  there  are 
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opportunities  in  methods  that  aid  intermodal  operators  in  determining  which  group 
of  routes  and  services  to  purchase  from  terminal  and  drayage  operators  on  behalf  of 
a  shipper  to  minimize  cost  and  delivery  time  frame. 

Zografos  and  Regan  discuss  challenges  in  intermodal  freight  transportation  in 
Europe  and  the  United  States  and  the  need  to  develop  information  and  communica¬ 
tion  technologies  to  aid  decision  making  in  a  deregulated,  free  enterprise  transporta¬ 
tion  market  [119].  They  reiterate  the  need  for  planning  methods  to  aid  3PL  firms 
in  purchasing  individual  services  from  separate  providers  that  can  act  as  integrated 
multimodal  transportation  options. 

The  “Multimodal  Corridor  and  Capacity  Analysis  Manual”  suggests  the  use 
of  multimodal  corridors  for  which  capacities  are  calculated  in  aggregate  [105].  One 
inadequacy  of  this  approach  is  that  any  change  or  disruption  in  a  link  within  the 
aggregated  corridor  capacity  estimate  would  require  a  separate  capacity  analysis. 

Other  network  capacity  estimation  research  efforts  lack  the  ability  to  model 
multimodal  freight  transport  capacity  and  allocation.  Asakura  estimates  the  capac¬ 
ity  of  urban  transportation  networks  given  user  route  preferences  using  a  bi-level 
programming  formulation  [8].  This  bi-level  programming  approach  is  also  used  by 
Park  and  Regan  to  estimate  the  capacity  of  multimodal  transportation  systems  al¬ 
though  the  model  is  not  solved  due  to  the  complexity  of  the  formulation  [96].  Yang 
et  al.  formulate  a  similar  network  capacity  problem  which  further  accounts  for 
the  potential  for  land  development  in  determining  route  selection.  Once  again  this 
approach  relies  on  a  bi-level  programming  formulation  for  which  a  heuristic  is  pro¬ 
posed  [115].  Morlok  and  Riddle  rely  upon  the  MAXCAP  problem  formulation  to 
estimate  transportation  system  capacity.  This  model  uses  aggregated  routes  based 
on  a  priori  information  of  commodity  flow  between  origin-destination  pairs  [88] .  Ge 
et  al.  use  a  bi-level  programming  formulation  to  determine  the  effect  of  reserve  ca¬ 
pacity  on  a  road  network  while  accounting  for  traveler  information  [49].  Methods 
not  relying  upon  bi-level  programming  to  plan  system  capacity  are  preferred  since 
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bi-lcvel  programming  problems  require  development  of  specialized  algorithms.  In  the 
case  of  mixed-integer  bi-level  programming  problems  solving  the  relaxation  does  not 
provide  a  valid  bound  on  the  optimal  solution  to  the  original  problem  causing  many 
difficulties  in  even  obtaining  feasible  solutions  to  capacity  problems  formulated  as 
bi-lcvel  programming  problems  [15,86]. 

The  preceding  paragraphs  provide  strong  evidence  that  simply  extending  cur¬ 
rent  modeling  techniques  is  not  sufficient  to  overcome  the  inherent  inadequacies  in 
current  capacity  modeling  approaches.  Therefore,  a  new  capacity  modeling  method¬ 
ology  is  required  in  order  to  improve  the  fidelity  and  accuracy  of  network  capacity 
estimates.  In  the  context  of  freight  transportation,  such  an  approach  must  capture 
important  system  parameters  like  single  trip  vehicle  capacity,  length  of  planning 
horizon  and  the  size  of  the  vehicle  fleet.  These  parameters,  and  others,  drive  net¬ 
work  capacity.  Allocation  of  this  capacity  should  provide  network  paths  efficiently 
using  fleet  capacity. 

Minimizing  the  data  collection  required  is  highly  desirable  since  data  collection 
is  often  difficult,  time-consuming,  and  expensive.  Reusing  data  that  are  already 
available  eliminates  the  need  to  collect  and  analyze  additional  data  sources  which 
may  support  no  other  purpose.  This  article  explicitly  addresses  all  of  these  factors 
and  develops  a  new  capacity  modeling  approach  for  capacity  planning.  Our  novel 
approach  is  then  extended  in  order  to  model  capacity  given  an  arbitrary  number  of 
distinct  transportation  modes. 

Given  a  fleet  of  vehicles  from  across  multiple  distinct  modes  of  transportation, 
single  source  and  sink,  and  planning  horizon  length,  consider  the  problem  of  deter¬ 
mining  the  multimodal  fleet  mix  and  path  which  will  provide  the  required  capacity 
at  minimum  cost  based  only  upon  operational-level  data  inputs. 

We  define  operational  data  as  those  data  that  are  available  without  requiring 
additional  data  collection.  Examples  would  be  engineering  and  operating  specihca- 
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tions  of  the  vehicle  fleets,  estimates  (in  tons)  of  the  freight  to  be  transported,  and 
the  length  of  the  time  window  for  transportation  planning. 

We  demonstrate  how  rough-cut  capacity  planning  (RCCP)  in  freight  networks 
is  modeled  using  a  modification  of  the  Resource  Constrained  Shortest  Path  Problem 
(RCSP).  The  basic  Shortest  Path  Problem  (SPP)  is  modified  by  adding  a  knapsack 
constraint  on  resource  consumption.  The  RCSP  was  originally  proposed  by  Handler 
and  Za.ng  in  1980  [54],  Define  the  following  sets  and  variables: 


N  =  {1,  2,  a  set  of  uniquely  labeled  nodes 

A  =  {( i,j )  :  i,  j  e  N}  a  set  of  directed  arcs  that  define 

adjacencies  for  the  nodes  in  N 
G  =  (TV,  A)  is  a  graph  defined  by  the  nodes  in 

N  and  the  arcs  in  A 

R  is  the  quantity  of  resource  available 

s  is  the  label  assigned  to  the  source 

node  where  supply  is  located 

t  is  the  label  assigned  to  the  destination 

node  where  demand  is  located 
Cij  is  the  cost  of  traversing  arc  (i,  j) 

Xij  =  1  if  arc  ( i,j )  is  included  in  the  path 

and  0  otherwise 


(3.1) 


The  RCSP  is  formulated  as: 
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Cij'X'ij 


(3.2a) 


Minimize 

(i,j)eA 

s.t.  xij-  Y  Xji  =  b(i)  Vi  e  N  (3.2b) 

j:(i,j)eA  j’-(j,i)eA 

CijXij  <  R  (3.2c) 

(i,j)eA 

b(s )  =  1, b(t)  =  —1, b{i)  =  0  V  i  ^  s}t  (3. 2d) 

In  this  formulation,  4.2a  is  the  objective  function  that  minimizes  cost,  4.2b 
are  the  flow  balance  constraints,  4.2c  are  the  knapsack  constraints  on  resources,  and 
4. 2d  are  the  flow  forcing  constraints. 

Applications  of  the  RCSP  to  transportation  include  solving  vehicle  routing 
problems  [43],  routing  a  homogeneous  vehicle  fleet  [11],  in  such  diverse  contexts 
as  wastewater  treatment  in  Kuwait  and  in  ensuring  compliance  with  government 
standards  for  thermal  resistance  in  building  design  [41]. 

The  RCSP  is  at  least  as  hard  as  NP-complete  problems  as  first  shown  by 
Handler  and  Zang  in  1980  [5,48,54,75].  RCSP’s  can  be  classified  according  to 
four  basic  criteria  [64]:  resource  accumulation,  path  structural  constraints,  objective 
function,  and  the  nature  of  the  underlying  network. 

Resources  may  be  accumulated  either  on  the  arcs  or  at  the  nodes.  Accumula¬ 
tion  of  resources  on  the  arcs  uses  resource  vectors  while  accumulation  at  the  nodes 
uses  resource  intervals.  These  two  methods  for  resource  accumulation  are  equivalent 
if  there  is  no  lower  bound  on  resource  consumption.  Path  structural  constraints  cap¬ 
ture  node- arc  adjacencies,  branch-and-price  rules,  and  admittance  of  any  non-simple 
paths  to  the  solution.  Objective  functions  are  classified  as  single-objective  or  multi¬ 
objective  depending  upon  the  number  of  objective  components  being  considered. 
Finally,  the  underlying  network  is  classified  as  cyclic  or  acyclic.  Networks  containing 
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cycles  have  an  infinite  number  of  paths  and  the  solutions  may  be  unbounded.  How¬ 
ever,  cyclic  networks  can  be  transformed  into  acyclic  networks  via  discretization  to 
prevent  unbounded  solutions  [64], 

We  now  turn  our  attention  from  classifying  RCSP  problems  to  classifying  RCSP 
solution  methods.  Garcia,  and  Mehlhorn  and  Ziegelmann  have  provided  classification 
schemes  for  RCSP’s  [47,82]  which  are  broadly  categorized  into  path  ranking,  labeling, 
and  Lagrangian  Relaxation  methods. 

Path  ranking  methods  are  based  upon  the  work  of  Eppstein,  and  of  Hoffman 
and  Paveley  who  first  studied  the  k-shortest  path  problem  [42,59,72],  Such  methods 
generally  operate  by  building  paths  and  then  ordering  them  according  to  some  cri¬ 
teria  such  as  cost  or  length.  In  RCSP’s  the  shortest  path  often  breaks  one  or  more 
of  the  knapsack  constraints.  Selecting  an  alternative  path  (from  the  rank  ordered  k- 
shortest  generated)  which  is  also  feasible  in  the  knapsack  constraints  is  one  approach 
using  k-shortest  paths. 

Labeling  methods  operate  by  searching  for  non-dominated  paths.  Each  path 
is  labeled  with  both  its  cost  and  resource  consumption.  When  an  alternate  path 
is  found  that  is  superior  to  the  incumbent  with  respect  to  either  cost  or  resource 
consumption  and  at  least  as  good  in  the  remaining  criterion,  then  the  incumbent  is 
replaced  by  the  alternate  path.  Irnich  and  Villeneuve  provide  a  labeling  method  for 
k-cyclc  elimination  for  k  >  3  [65]. 

Lagrangian  Relaxation-based  solution  methods  solve  a  relaxation  of  the  original 
problem  by  relaxing  the  knapsack  constraint  associated  with  the  resource  consump¬ 
tion,  and  then  working  to  close  the  duality  gap  between  the  relaxed  linear  program¬ 
ming  primal  and  dual  solutions.  Mehlhorn  and  Ziegelmann  have  created  a  RCSP  soft¬ 
ware  package  called  “Constrained  Network  Optimization  Package”  which  operates 
based  on  Lagrangian  Relaxation  and  closing  the  primal-dual  gap  [11,18,21,100,113]. 
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Heuristics  other  than  Lagrangian  Relaxation  have  also  been  successfully  ap¬ 
plied  to  the  RCSP  and  similar  problems.  Particle  Swarm  Optimization  was  applied 
to  the  Resource  Constrained  Project  Scheduling  Problem  by  Zhang  et.  al.  [116].  Has- 
sin  developed  a  fully  polynomial  approximation  bounding  scheme  for  the  RCSP  that 
was  later  generalized  by  Lorenz  and  Raz  to  accommodate  both  positive  and  negative 
arc  weights  [57,78].  Hu  et  al.  develop  an  improved  version  of  the  Ant  System  (AS) 
originally  proposed  by  Dorigo  for  solving  the  RCSP  [63].  The  proposed  AS  relies 
upon  modified  heuristic  information  and  pheromone  updates  to  drive  convergence 
to  paths  meeting  the  resource  constraints  [33] . 

Zhu  and  Wilhelm  developed  a  three  step  approximation  algorithm  that  is  ap¬ 
plied  iteratively  to  obtain  successively  tighter  bounds  on  the  optimal  solution  [117]. 
Step  one  preprocesses  the  graph  to  remove  nodes  and  arcs  which  cannot  be  contained 
on  any  feasible  path.  Step  two  expands  the  reduced  graph  based  upon  resource  avail¬ 
ability  and  step  three  solves  the  unconstrained  shortest  path  problem.  Arc  weights 
are  adjusted  and  the  algorithm  is  repeated  until  sufficiently  tight  bounds  are  deter¬ 
mined.  There  are  many  other  research  efforts  that  have  been  devoted  to  developing 
fast  running  heuristics  for  the  RCSP  [10,39,57,75]. 

Jepsen  et  al.  provided  the  first  solution  method  for  the  ERSCP  based  on 
branch- and-cut.  In  order  to  strengthen  the  formulation  a  new  class  of  valid  inequality 
called  generalized  capacity  inequalities ,  is  developed.  This  new  class  of  inequality  is 
combined  with  knapsack  constraints  and  subtour  elimination  constraints  which  are 
all  introduced  as  cutting  planes  in  the  branch-and-cut  algorithm  [67]. 

No  source  reviewed  in  the  literature  provides  a  time  scalable  method  of  multi¬ 
modal  freight  transportation  network  capacity  estimation.  Scalability  refers  to  the 
ability  of  the  method  to  adapt  capacity  estimates  to  an  arbitrarily  defined  planning 
horizon.  Furthermore,  the  limited  number  of  capacity  estimation  methods  are  based 
on  tactical  level  considerations  and  are  not  well  suited  to  estimating  overall  system 
capacity  without  tactical  level  data  which  may  not  be  available  [96]. 
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3.1.1  Assumptions.  For  the  purposes  of  this  research  we  assume  that 
capacity  is  expressed  in  terms  of  tons.  This  assumption  allows  us  to  use  million-ton 
miles  (MTM)  as  the  resource  for  the  knapsack  constraints  as  we  demonstrate  in  the 
next  section. 

We  do  not  consider  weight  and  balance,  dimensional  packing,  or  individual 
vehicle  routing  as  these  problems  are  well  represented  in  the  literature  [3,13,20,70,80]. 

Finally,  we  assume  that  scheduling  is  completed  in  near-term  planning.  RCCP 
is  targeted  at  mid-range  planning  and  answering  the  question  “Does  the  system  have 
sufficient  capacity  to  meet  the  demand.”  Once  demand  and  capacity  are  equilibrated 
through  RCCP,  scheduling  algorithms  and  heuristics  are  used  to  build  final  delivery 
schedules  [58,97]. 

3.2  Modeling  Approach 

The  modeling  approach  established  in  this  section  treats  each  transportation 
mode  as  a  separate  layer  of  a  network.  Each  network  layer  shares  a  common  node 
set  with  all  other  network  layers.  The  layers  differ  by  the  arc  sets  which  define  node 
adjacencies  within  each  layer.  In  order  to  ensure  that  the  available  capacity  within 
each  mode  of  transportation  is  respected  by  the  modeling  approach,  the  RCSP  is 
used  to  capture  the  capacity  constraints  which  are  defined  for  each  layer.  These 
capacity  constraints  will  be  defined  in  the  following  sections. 

The  RCSP  is  preferred  to  the  Multicommodity  Network  Flow  Problem  (MCNF) 
presented  by  Ahuja  et  al.  since  the  MCNF  uses  bundling  constraints  on  total  arc 
flow  shared  by  all  commodities  and  cannot  capture  individual  capacity  constraints 
for  the  various  modes  of  transportation  [2],  By  contrast,  the  knapsack  constraints 
on  resources  can  be  used  to  constrain  resource  usage  through  summing  arc  weights 
along  a  network  path. 

Modeling  the  rough-cut  capacity  problem  for  multimodal  transportation  re¬ 
quires  the  development  of  some  model  parameters.  First,  we  formulate  the  single 
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mode  case  as  a  RCSP  problem  and  then  generalize  the  result  to  an  arbitrary  number 
of  distinct  transportation  modes.  In  a  physical  network,  each  pairwise  adjacent  set 
of  nodes  has  associated  with  it,  a  distance  label.  The  graph  is  not  required  to  be 
complete.  That  is,  not  all  nodes  in  the  network  are  pairwise  adjacent.  For  example, 
consider  a  network  of  four  shipping  and  receiving  locations.  Let  node  1  be  a  source 
node,  node  4  be  a  sink  node,  and  nodes  2  and  3  be  transshipment  nodes.  The  arc 
labels  represent  the  Euclidean  distance  between  nodes.  Such  a  network  is  shown  in 
Figure  3.1: 


Network  capacity  can  be  viewed  as  the  ability  of  a  network  to  provide  a  re¬ 
quired  level  of  throughput  within  a  given  time  window.  Since  network  flow  may  be 
continuous,  the  amount  of  throughput  must  be  defined  for  a  specified  time  period 
in  order  to  provide  a  basis  for  measurement  and  comparison.  Capacity  estimates 
in  freight  transportation  networks  are  conveniently  expressed  in  units  of  million  ton 
miles  per  day  (MTM/D).  This  standard  measure  of  vehicle  fleet  performance  is  the 
standard  for  air  mobility  but  may  be  calculated  for  any  fleet  of  vehicles  [91].  Implicit 
in  this  definition  of  capacity  is  that  the  “cost”  associated  with  adding  capacity  to  a 
route,  concerns  distance,  time,  and  required  capacity  (expressed  in  millions  of  tons). 

In  figure  3.1,  for  example,  the  distance  between  node  1  and  node  2,  d12  = 
d'2i  =  8.  If  five  million  tons  of  freight  must  be  transported  between  nodes  1  and  2 
at  a  cost  of  $2  per  MTM  then  we  can  re-weight  the  arcs  in  the  network,  expressing 
weights  as  the  cost  for  the  required  million-ton-miles  (MTM)  allocated  to  an  arc 


by  multiplying  each  dij  by  the  number  of  tons  of  freight  (in  millions)  required  to 
be  moved,  m,  and  the  cost  per  MTM  h.  This  label  must  then  be  weighted  by  the 
estimated  cost  per  MTM  for  the  mode  of  transportation.  The  final  arc  weight  is 
then  given  by:  csij  =  dtl  *m*  h.  The  arc  weights  now  express  a  composite  quantity 
which  addresses  cost,  distance,  and  required  capacity.  In  order  to  assign  capacity  to 
arc  Cij/h  MTM  are  required.  The  reweighted  network  is  labeled  as  shown  in 

Figure  3.2: 


The  number  of  MTM  available,  r,  is  a  function  of  the  transportation  fleet 
characteristics  and  the  planning  horizon  length  (T).  This  can  be  expressed  as:  r  = 
V  *  MT  *  T  where  MT  is  the  million-ton-mile  per  day  capability  of  a  single  vehicle 
and  V  is  the  number  of  vehicles  in  the  fleet.  By  inspection  the  shortest  path  through 
the  example  network  is  1  — >  2  — y  4.  This  path  requires  45  MTM  at  a  cost  of  $90.  If 
our  fleet  cannot  provide  the  required  MTM,  then  the  problem  is  infeasible  unless  we 
add  assets  to  the  fleet  (increase  V )  or  add  days  to  the  extend  the  planning  horizon 
length  (T). 

The  parameter  MT  is  calculated  based  on  the  following  operational  level  data: 
“UTE”  rate,  payload,  block  speed,  and  productivity.  UTE  rate  is  the  average  his¬ 
torical  “surge”  capability  of  a  single  vehicle  in  hours  per  day.  This  represents  the 
maximum  availability  of  each  vehicle  per  day.  Down  time  clue  to  scheduled  main¬ 
tenance  and  reliability  issues  are  assumed  to  be  captured  in  this  quantity.  Vehicle 
payload  is  the  single  trip  capacity  of  a  vehicle  in  tons.  Vehicle  block  speed  is  the 
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average  door-to-door  vehicle  speed.  Productivity,  refers  to  the  average  historical 
percentage  of  non-empty  trips  over  all  vehicles.  This  parameter  accounts  for  repo¬ 
sitioning  of  vehicles  within  the  fleet  when  they  are  moved  while  empty  in  order  to 
pick  up  a  load  of  freight.  Now  MT  is  calculated  as  MT  =  MTM/D  per  asset  = 
(f/TPNPayload*Block  Speed*Productivity)/l,  000, 000.  The  single  mode  formula¬ 
tion  can  now  be  presented  building  upon  the  development  of  these  parameters. 

3.2.1  Single  Mode  Formulation.  To  simplify  the  presentation  of  the  single 
mode  resource  constrained  shortest  path  formulation  the  following  parameters  and 
vectors  are  defined:  N  =  {1,2,  is  the  set  of  nodes  representing  shipping  loca¬ 

tions,  C  =  {( l,v )  :  (l,v)  is  an  arc  in  the  graph}  is  the  set  of  feasible  directed  arcs 
defined  over  the  elements  in  N,  G  =  (N,  C )  is  the  graph  consisting  of  nodes  in  N 
and  arcs  in  C,  and  E  is  the  n  x  \C\  node-arc  adjacency  matrix  where  each  column  of 
E  has  a  1  in  the  Ith  row  and  a  —1  in  the  Vth  row  corresponding  to  adjacency  matrix 
element  Eiv.  Now  define  Eiv  along  with  the  following  sets,  matrices,  and  scalars  as: 
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Eij  - 


1  if  l  G  N  is  the  tail  of  arc  j  G  C 

—  1  if  Z  G  N  is  the  head  of  arc  j  e  C  V  l  G  iV,  j  G  C 

0  otherwise 


1  if  ESj  =  1 
0  otherwise 


V  iec 


g-f(j)  = 


1  if  Etj  =  -1 

0  otherwise 


V  jeC 


Q  = 


(3.3) 


h 


the  cost  (in  dollars)  per  MTM 


r  : 


the  number  of  MTM  available 


B_s  :  the  1  x  \C\  vector  of  arc  weights 


x  :  the  1  x  \C\  vector  of  decision  variables 


J  = 


B_s 

~h~ 


The  single  mode  formulation  can  now  be  expressed  as: 
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Minimized  s  x 


(3.4a) 

subject  to :Ex  =  0  (3.4b) 

Qx  =  1  (3.4c) 

Jx  <  r  (3.4d) 

x  E  B|c|  (3.4e) 


By  comparison  with  the  generic  RCSP  formulation  presented  in  section  3.1, 
this  formulation  is  also  a  RCSP  problem  where  5.2a  is  the  objective  function,  5.2b 
are  the  transshipment  portion  of  the  flow  balance  constraints,  5.2c  are  the  source 
and  sink  portion  of  the  flow  balance  constraints,  5. 2d  are  the  knapsack  constraints, 
and  5.2e  are  the  flow  forcing  constraints. 

3.2.2  Multimodal  Formulation.  Generalizing  the  single  mode  model  to 
incorporate  multiple  transportation  modes  requires  introduction  of  several  other  sets 
and  indices  to  aid  in  describing  not  only  the  intra-mode  arcs  for  each  transportation 
mode  but  also  the  inter-mode  arcs.  To  develop  the  generalization,  an  identical  “copy” 
of  the  physical  node  set  is  created  for  each  mode  of  transportation.  Although  the 
nodes  for  each  mode  of  transportation  are  identical,  the  edge  set  for  each  mode  of 
transportation  may  be  distinct  in  order  to  represent  different  node  adjacencies  that 
might  exist  depending  upon  the  mode  of  transportation.  Furthermore,  some  nodes 
may  act  as  intermodal  transshipment  points  at  which  freight  may  be  transfered 
between  modes  of  transportation.  A  port  which  has  both  a  rail  head  and  a  sea  port 
may  be  one  such  combination.  These  intermodal  transfers  require  definition  of  new 
sets  of  intermodal  arcs.  An  example  is  provided  following  the  development  of  the 
multimodal  RCSP  formulation. 
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We  first  define  the  following  sets  and  parameters  to  generalize  the  single  mode 
for  addressing  multimodal  freight  transportation.  Let  P  —  {1,  2, . . .  ,p}  be  the  set 
of  available  transportation  modes  and  k  G  P  be  be  the  index  for  these  modes.  Also 
N  =  {1,2 , rn}  be  the  set  of  nodes,  as  in  the  single-mode  formulation  with  i  G  N 
being  the  index  for  the  nodes.  Expanding  the  original  graph  as  described  in  the 
previous  paragraph  results  in  a  node  set  with  rip  nodes.  The  indices  for  the  node  set 
of  the  expanded  graph  are  /,  v  G  {1, . . . ,  np}.  Let  Nk  =  {l+n(k—  1), . . . ,  n+n(k— 1)} 
be  the  node  set  for  mode  k  in  the  expanded  graph.  Finally,  let  Ml  =  {i, . . .  ,i  +  n(k  — 
1), . . . ,  i  +  n(p  —  1)}  be  the  node  set  in  the  expanded  graph  corresponding  to  node 
i  in  the  original  graph.  Define  the  additional  sets,  and  matrices  used  in  formulating 
the  multimodal  RCSP: 


Ck  =  {(l,v)  :  (l,v)  is  an  intra-mode  arc  in  the  graph  for  mode  k } 

Dl  =  {(l,v)  :  (l,v)  is  an  inter-mode  arc  in  the  graph  for  node  i} 

V  =  U^':  n°de  set  for  the  expanded  graph 

k 

A  =  ^|J  Ck^j  |J  D^j  :  the  arc  set  for  the  expanded  graph  (3-5) 

Gk  =  ( Nk,Ck ):  the  graph  for  mode  k 
Hl  =  (M\  Dl):  the  graph  for  node  i 
T  =  (p,  A):  the  expanded  graph 

Now  we  define  additional  indices.  Let  u  G  {1, ...,  |A|}  be  the  index  for  the  arc 
set  of  the  expanded  graph,  mk  G  {1,  ...|Cfc|}  be  the  index  of  the  intra-mode  arc  set 
for  for  each  transportation  mode  k,  and  zl  G  {1, ...,  \Dl\}  be  the  index  of  the  inter¬ 
mode  arc  set  for  each  node  i.  Also  let  s  and  t  be  the  indices  of  the  source  and  sink 
nodes  in  the  original,  non-expanded  graph  respectively.  With  these  sets  and  indices 
we  can  now  develop  the  multimodal  formulation  using  the  matrices  defined  below. 
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1 

if  /  G  Nk  is  tail  of  arc  mk  e  Ck 

rpk  _ 

< 

-1 

if  /  G  Nk  is  head  of  arc  mk  G  C 

0 

otherwise 

1 

if  /  G  Ml  is  tail  of  arc  zl  G  Dl 

F*  — 

r  l,zx 

< 

-1 

if  /  G  Ml  is  head  of  arc  zl  G  Dl 

0 

V 

otherwise 

E  = 

[E1:.E2: . 

. . : E'p } 

F  = 

[EV:E2:. 

. .  \Fn] 

A  = 

[e:f] 

qs(u)  =  Y  Ai,uVu 

ieMs 
Ai,u  =  1 

q-t(u)  =  Y  A  i,uVu 

leM* 

Al,u=~l 


Q  = 


q-t 


(3.6) 


Additionally,  we  define  the  following  vectors  and  matrices.  Note  that  in  the 
definitions  below,  the  notation  B'  is  taken  to  be  the  matrix  transpose  of  B. 
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h  :  the  1  x  p  matrix  of  mode  transportation  costs  in  dollars  per  MTM 

r  :  the  1  x  p  vector  of  MTM  available  in  mode  k 

B_sk  :  the  1  x  \Ck\  vector  of  arc  weights  for  mode  k 


BP  :  the  1  x  \D'l\  vector  of  arc  weights  for  node  i 

B_s  =  [B  s^-.B  s2\ . . .  \B_sp] 

B_t  =  [B_tr:B_t2\  .  .  .  \B_tn } 

x  :  the  1  x  |A|  vector  of  decision  variables 

\B  s'] 


B'  = 

B_t' 


J  =  [JV.  0] 


Using  the  preceding  developments  we  can  now  formulate  the  multimodal  RCSP 


as: 


Minimize  Bx 
subject  to:  Ax  =  0 

Qx  =  1  (3.8) 

Jx  <  r 
x  e  b|A| 

3.2.2. 1  Dimensional  Analysis  &  Problem  Structure.  In  many  real- 
world  math  programming  applications,  the  constraint  matrix  exhibits  special  struc¬ 
ture  which  may  be  exploited  to  quickly  solve  large  problems  [77].  In  fact,  many 
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commercial  solvers  perform  preprocessing  analysis  of  problems  in  order  to  determine 
if  such  structure  is  present  in  a  given  problem  instance.  These  solvers  will  exploit  the 
problem  structure  internally  without  intervention  by  the  user.  The  constraint  matrix 
generated  in  the  previous  section  exhibits  a  bordered  angular  structure  exploitable 
through  cross  decomposition  which  is  a  combination  of  Dantzig- Wolfe  and  Benders’ 
decomposition.  Before  presenting  the  entire  constraint  matrix,  we  begin  by  showing 
consistency  of  dimensions  among  all  constraint  matrix  subcomponents. 

Each  matrix,  has  a  column  for  each  arc  contained  in  the  intra-mode  graph 
of  mode  k  .  Each  Fl  matrix  has  a  column  for  every  arc  in  the  inter-mode  graph  on 
node  i.  Both  Ek  and  Fk  have  a  row  for  every  node  in  the  expanded  graph  as  shown 
in  the  figures  below.  Note  that  Ek  will  have  rows  of  zeros  for  all  nodes  l  not  in  mode 
k  \/k  G  P.  The  structure  of  each  matrix  Ek  results  in  the  matrix  E  having  specal 
structure  as  discussed  in  the  proceeding  paragraph. 


The  node-arc  incidence  matrix,  A,  for  the  expanded  graph  has  a  row  for  each 
node  and  a  column  for  each  inter-modal  and  intra-modal  arc.  A  comprises  both  E, 
the  matrix  formed  by  horizontally  concatenating  all  matrices  Ek  and  F,  the  matrix 
formed  by  horizontally  concatenating  all  matrices  Fl. 


The  following  graphic  shows  the  special  block  diagonal  structure  of  the  E  ma¬ 
trix.  This  structure  is  due  to  the  way  the  expanded  multimodal  graph  is  constructed 
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by  creating  a  “copy”  of  the  nodes  with  unique  labels  ofr  each  transportation  mode. 
Each  E  in  the  structural  representation  of  E  below  depicts  a  block  of  non-zero  ele¬ 
ments  and  the  0  entries  depict  appropriately  dimensioned  blocks  of  0  entries. 


EC* 


k 


E 

0 

0 

0 

0 

E 

0 

0 

0 

0 

E 

0 

0 

0 

0 

E 

The  vectors  in  the  constraint  matrix  that  ensure  exactly  one  source  and  one 
sink  node  are  chosen  have  a  single  row  and  a  have  column  for  every  arc  in  the 
expanded  graph. 


qs 


1 


|A[ 

qt 


l 


|A| 


Q 


2 


The  dimensional  analysis  for  the  B_s  and  B_t  is  identical  to  the  analysis  for  E 
and  F.  These  vectors  are  provided  below  for  completeness  along  with  the  vector, x, 
of  decision  variables. 


The  vector,  d,  has  a  column  for  each  arc  in  the  expanded  graph.  The  matrix, 
./]  has  a  row  for  each  intra-modal  arc  in  the  expanded  graph  and  a  column  for  each 
mode  of  transportation  appending  an  appropriately  dimensioned  matrix  of  zeros  to 
J\  yields  the  matrix  J  with  a  column  for  each  arc  in  the  expanded  graph  and  a  row 
for  each  mode  of  transportation. 
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The  graphic  below  depicts  the  structure  of  the  complete  constraint  matrix.  The 
dimensions  of  the  complete  constraint  matrix  are  (p(n  +  1)  +  2)  x  |A|  are  obtained 
by  summing  the  component  rows  and  columns. 


This  modeling  approach  represents  a  capacity  planning  at  a  Rough-Cut  level 
of  fidelity  since  it  requires  only  operational  level  inputs.  This  method  of  capacity 
estimation  is  analogous  to  the  “Capacity  Planning  based  on  Overall  Factors”  which 
uses  historical  data  to  perform  mid-range  to  long-range  planning  in  the  context  of 
a  hierarchical  production  planning  process.  Readers  unfamiliar  with  Hierarchical 
Production  Planning  are  encouraged  to  reference  Hax  and  Candea  or  Bitran  and 
Hax  [19,58]. 

One  feature  of  the  model  that  is  important  to  note  is  that  it  will  always  return 
the  least  cost  feasible  shipping  option  for  the  required  capacity.  This  means  that  if 
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a  single  mode  option  is  least  expensive  and  the  fleet  size  for  the  mode  is  capable  of 
providing  the  required  number  of  million-ton-miles,  then  the  solution  returned  will 
be  a  single  mode  solution. 

3.3  Example  Problem.  Formulation 

In  this  section  we  explore  the  nature  of  the  methodology  presented  in  the 
previous  section  by  formulating  a  small  example  problem  to  aid  in  methodology 
conceptualization  and  to  show  empirically,  the  utility  of  the  modeling  approach  de¬ 
veloped.  An  interesting  feature  of  the  modeling  approach  is  that  it  accounts  for  the 
intermodal  transfer  costs.  Due  to  this  fact,  the  model  will  perform  mode  selection  if 
a  single  mode  of  transportation  is  both  feasible  and  least  expensive.  Alternatively, 
if  an  integrated  multimodal  route  is  less  expensive,  or  if  it  is  the  only  feasible  option 
due  to  a  lack  of  capacity  in  other  transportation  modes,  then  the  model  will  select 
an  integrated  multimodal  route. 

The  formulation  presented  below  was  constructed  to  demonstrate  a  situation 
where  an  integrated  multimodal  route  solution  is  optimal.  This  same  problem  is 
modified  slightly  (the  intermodal  transfer  weights  are  increased)  in  order  to  demon¬ 
strate  the  use  of  the  model  for  single  mode  selection.  The  final  characteristics  demon¬ 
strated  using  this  small  example  problem  are  the  affect  on  capacity  (and  the  resulting 
optimal  solution)  due  to  changes  in  capacity  and  planning  horizon  length. 

3.3.1  An  Example  Formulation.  Beginning  with  the  single-mode  graph 
presented  in  section  3.2,  we  expand  the  graph  to  include  two  additional  modes  of 
transportation.  Each  of  the  transportation  modes  has  a  different  node  adjacency 
structure.  The  intra-modal  graphs  of  the  two  additional  transportation  modes  are 
shown  below.  These  graphs  represent  G'1,  G2,  and  G3. 

Now  in  order  to  define  the  complete  expanded  graph,  the  nodes  must  be  rela¬ 
beled  according  to  v  =  i  +  n(k—  1)  and  the  inter-modal  graphs  defined  on  each  node 
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in  the  original  graph  must  be  generated.  The  inter-modal  graphs  are  shown  below. 
These  graphs  shown  in  Figure  3.5  (from  left  to  right)  are  the  graphs  H1,  H 2,  H 3, 
and  HA\ 


By  combining  the  intra-mode  and  inter-mode  graphs  it  is  now  possible  to 
specify  a  single  graph  for  use  in  formulating  the  multimodal  RCSP  as  shown  in 
Figure  3.6. 

The  sets  used  to  formulate  the  preceding  graph  are  provided  below  for  com¬ 
parison  and  reference.  In  this  example  n  —  4,  p  —  3. 
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i  6  {1,2, 3, 4},  *£{1,2,3} 

l,v  £  {1,2,3,4,5,6,7,8,9,10,11,12} 


s  —  1,  t  —  4 

iV1  =  {1,2,3,4},  IV2  =  {5,6,  7,  8},  N3  =  {9, 10, 11, 12} 

M1  =  {1,  5,  9},  M2  =  {2,6, 10},  M3  =  {3,7, 11},  M4  =  {4,8, 12} 
C1  =  {(1,  2) ,  (2, 1) ,  (2,  3) ,  (3, 2) ,  (2, 4) ,  (4,  2)} 

C2  =  {(5,7),  (7, 5),  (6,  7),  (7, 6)} 

C3  =  {(9,11),  (11,  9),  (10, 11),  (11, 10)} 

D1  =  {(  1,9),  (9,1)} 

D2  =  {(  2,  6),  (6, 2)},  D3  =  {(3,7),  (7,  3),  (7, 11),  (11,  7)},  D4  =  {} 


77  =  {1,2,  3, 4,  5,  6,  7, 8,  9, 10, 11, 12} 

'  (1,  2) ,  (2, 1) ,  (2,  3) ,  (3,  2) ,  (2, 4) ,  (4,  2) ,  (5,  7) ,  (7,  5) , 
A  =  (6, 7) ,  (7,  6) ,  (9, 11) ,  (11,  9) ,  (10, 11) ,  (11, 10) ,  (1,  9) , 

(9,1),  (2,  6),  (6, 2),  (3,  7),  (7,  3),  (7, 11),  (11,  7) 


(3.9) 


Indices  relying  upon  the  above  set  definitions  are  u,  mk ,  and  zl.  These  indices 
and  the  matrices  they  index  are: 


rl-10  o  0  O-i 
-11  1-10-1 
0  0-1100 
o  o  o  0-11 
oooooo 
oooooo 
oooooo 
oooooo 
oooooo 
oooooo 
oooooo 
L000000J 
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e2 


r  o  o  o  o-i 
oooo 
oooo 
oooo 
1-10  0 
0  0  1-1 
-11-11 
oooo 
oooo 
oooo 
oooo 
L  0  0  0  0  J 


E3 


r  o  o  o  o  -i 
oooo 
oooo 
oooo 
oooo 
oooo 
oooo 
oooo 
1-10  0 
0  0  1-1 
-11-11 
L  0  0  0  0  J 


E 


r  1  -1  0  0  0  0 

-11  1-10-1 
0  0-11  0  0 

0000-11 
0  0  0  0  0  0 

0  0  0  0  0  0 

0  0  0  0  0  0 

0  0  0  0  0  0 

0  0  0  0  0  0 

0  0  0  0  0  0 

0  0  0  0  0  0 

L  0  0  0  0  0  0 


0  0  0  0  0  0  0 

0  0  0  0  0  0  0 

0  0  0  0  0  0  0 

0  0  0  0  0  0  0 

1  -1  0  0  0  0  0 

0  0  1  -1  0  0  0 

-11-11  0  0  0 
0  0  0  0  0  0  0 

0  0  0  0  1  -1  0 

0  0  0  0  0  0  1 

0  0  0  0-11-1 
0  0  0  0  0  0  0 


0  1 

0 

0 

0 

0 

0 

0 

0 

0 

-1 

1 

0  j 


Fi  = 


r  1  -1! 
0  0 
0  0 
0  0 
0  0 
0  0 
0  0 
0  0 
-1  1 
0  0 
0  0 
0  0  J 


F2 


r  o  o-i 

1  -1 
0  0 
0  0 
0  0 
-1  1 
0  0 
0  0 
0  0 
0  0 
0  0 
0  0  J 


Fs 


r  o  o  o  o-i 
oooo 
1-10  0 
oooo 
oooo 
oooo 
-11  1-1 
oooo 
oooo 
oooo 
0  0-11 
L  0  0  0  0  J 
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Omitted  since  it  is  a  matrix  of  zeros 


F4  = 


F  = 


r  1  -1  00000 
-11  1-10-10 
00-11000 
0000-110 
0  0  0  0  0  0  1 

0  0  0  0  0  0  0 

o  o  o  o  o  o  -: 

0  0  0  0  0  0  0 

0  0  0  0  0  0  0 

0  0  0  0  0  0  0 

0  0  0  0  0  0  0 

Looooooo 


1-1000000' 
0  0  1  -1  0  0  0  0 

00001  -1  00 
00000000 
00000000 
00-110000 
0  0  0  0  -1  1  1  -1 

00000000 
-1  1000000 
00000000 
000000-11 
00000000. 


00000001 
00000000 
00000000 
00000000 
-10000000 
01  -1  00000 
1-1100000 
00000000 
0  0  0  1  -1  0  0  -1 
000001-10 
0  0  0-11-11  0 
00000000 


-1  000000-1 
01-10000 
0001-100 
0  0  0  0  0  0  0 

0  0  0  0  0  0  0 

0  -1  1  0  0  0  0 

0  0  0  -1  1  1  -1 

0  0  0  0  0  0  0 

1  0  0  0  0  0  0 

0  0  0  0  0  0  0 

0  0  0  0  0  -1  1 

oooooooJ 


qs  =  1000001000100010000000 


Q-t  = 

Q  = 


0000100000000000000000 

1000001000100010000000 

0000100000000000000000 


Before  the  formulation  can  be  completed,  the  distance  labels  applied  to  the 
multimodal  network  must  be  modified  by  accounting  for  the  parameters  associated 
with  time  and  shipping/transfer  costs.  Let  h  =  [0.4,  2,2]  be  the  vector  of  shipping 
costs  per  MTM.  There  is  an  entry  in  the  vector  for  each  of  the  three  transportation 
modes.  Similarly,  assume  the  remaining  parameters  are  specified  as  follows:  number 
of  vehicles  in  each  mode  is  given  by  the  vector  [30,21,20],  block  speed  in  miles 
per  hour  is  specified  by  [300,400,500],  UTE  rate  in  hours  per  day  by  [10,12,15], 
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single  trip  payload  for  each  vehicle  type  in  units  of  tons  is  given  by  the  vector 
[20,20,20],  and  the  productivity  as  an  average  historical  percentage  of  trips  where 
each  vehicle  type  travels  with  a  full  freight  load  as  [.5,  .5,  .5].  Additionally,  let  us 
assume  we  have  T  =  100  days  to  complete  the  freight  movement  and  assume  we 
wish  to  transport  ten  million  tons  of  freight  (j)  =  10.  Finally,  assume  the  intermodal 
transshipment  costs  are  $1  per  million  tons.  The  values  in  B_t,  the  inter-modal  arc 
weights,  are  calculated  based  upon  the  transfer  cost  per  MTM  and  the  number  of 
tons  (in  millions)  of  demand  (0).  An  intermodal  arc  from  mode  1  to  node  9  (SJ(1)), 
for  example  is  reweighted  as  B_t(  1)  =  h(  1)  *  (f>  =  1  *  10  =  10.  The  reweighed  graph 
represents  the  graph  for  ijj. 


Modifying  the  arc  weights  in  order  to  represent  capacity  in  terms  of  MTM , 
we  can  now  generate  the  B_s,  B_t,  d,  J\  and  J  matrices  and  and  the  r  vector. 
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r  = 


90  96  150 


B.s 


32  32  4  4  4  4  20  20  10  10  20  20  20  20 


BJt 


10  10  10  10  10  10  10  10 


d  = 


32  32  4  4  4  4  20  20  10  10  20  20  20  20  10  10  10  10  10  10  10  10 


Jl  = 


80  80  10  10  10  10  0  0  0  0  0  0  0  0 
000000  10  10  550000 

0  0  0  0  0  0  0  0  0  0  10  10  10  10  J 


80  80  10  10  10  10  0  0000  0  0  000000000 

0  0  0  0  0  0  10  10  5  50  0  0  000000000 

L  0  0  0  0  0  0  0  0  0  0  10  10  10  10  0  0  0  0  0  0  0  0  J 


This  completes  the  formulation  of  the  constraint  matrix  components,  ft  is  a 
trivial  exercise  to  verify  consistency  among  the  dimensions  and  formulate  the  entire 
constraint  matrix  and,  thus,  these  final  steps  are  omitted. 

3.3.2  Mode  Selection  &  Integration.  One  feature  of  the  modeling  approach 
important  to  note  is  that  it  will  always  return  the  least  cost  feasible  shipping  option 
for  the  required  capacity.  This  means  that  if  a  single  mode  option  is  least  expensive 
and  the  fleet  size  for  the  mode  is  capable  of  providing  the  required  number  of  million- 
ton-miles,  then  the  solution  returned  will  be  a  single  mode  solution. 

To  demonstrate  this  capability,  consider  the  example  formulation  of  the  pre¬ 
vious  section.  The  solution  is  given  by  the  single  modal  path  1  — >•  2  — >•  4.  If  modal 
transportation  costs  change  then  optimality  may  also  change.  Similarly,  if  resource 
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constraints  change  then  feasibility  may  change.  Either  of  these  changes  has  the 
potential  to  affect  the  optimal  solution. 

Now  assume  that  the  transportation  cost  vector  changes  to  h  =  [5,2,2],  The 
updated  values  for  Bt  are  shown  below  and  the  new  optimal  solution  is  given  and 
interpreted. 


B_s 


400  400  50  50  50  50  20  20  10  10  20  20  20  20 


The  optimal  solution  now  changes  to  the  bi-modal  path  5— >  7  — >  6  — >  2  — >4. 
This  path  is  interpreted  in  the  physical  network  as  beginning  at  node  1  in  mode  2 
and  adding  capacity  to  arc  (1,  3).  At  node  3  in  mode  2,  capacity  is  added  to  (3,  2)  in 
mode  2.  Capacity  is  transfered  from  node  2  in  mode  2  to  node  2  in  mode  1.  Finally, 
capacity  is  added  to  arc  (2,4)  in  mode  1. 

Consider  another  mode  price  vector  of  h  =  [5,  2,  0.5].  This  transportation  price 
structure  admits  a  truly  multimodal  optimal  solution  in  which  all  three  modes  of 
transportation  are  used  in  an  integrated  way  to  provide  the  required  capacity.  The 
updated  Bt  vector  for  the  new  price  structure  is: 


B_s 


400  400  50  50  50  50  20  20  10  10  5  5  5  5 


The  multimodal  solution  corresponding  to  this  value  of  h  is  9  — >■  11  — >■  7  — >■ 
6  — >  2  — >■  4.  Interpretation  of  this  multimodal  path  is  similar  to  that  for  the  bi- 
modal  path  and  is  left  to  the  reader  as  a  simple  exercise.  The  next  section  discusses 
computation  results  of  the  methodology. 


3-4  Results 

Problems  in  this  section  are  randomly  generated  according  to  several  param¬ 
eter  inputs:  Number  of  nodes,  number  of  modes,  density  of  arcs,  and  time  and 
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transportation  asset  availability.  The  problems  are  generated  using  MATLAB  and 
output  as  Mathematical  Programming  System  (MPS)  files.  The  MPS  file  format  is  a 
standard  modeling  format  which  is  recognized  by  major  math  programming  solvers 
like  LINGO  and  CPLEX.  Models  are  solved  using  LINGO  11.0  running  on  a  Toshiba 
Satellite  E105  laptop  computer  with  a  2.26  GHz  Intel  Core2  Duo  CPLT  and  4  GB  of 
RAM  running  64-bit  Microsoft  Windows  XP. 

3-4-1  Parameter  Selection  &  MPS  File  Generation  Issues.  Generating  the 
MPS  files  was  time  consuming.  In  fact,  MPS  file  generation  took  longer  that  solving 
the  problem  represented  by  the  MPS  file.  In  addition,  file  sizes  for  the  generated 
files  were  very  large.  The  largest  file  generated  contained  250  nodes  and  3  modes 
of  transportation  for  a  total  of  750  nodes.  The  MPS  file  size  for  this  problem  was 
approximately  28  MB.  The  authors  feel  that  the  generated  problems  are  of  reasonable 
size  to  demonstrate  the  validity  of  the  RCSP  formulation. 

3-4-2  Table  of  Results  for  Successfully  Generated  Problems.  Problems  were 
randomly  generated  using  90%  arc  density  and  three  modes  of  transportation.  The 
smallest  problem  has  150  total  nodes  (50  nodes  and  three  transportation  modes)  and 
probelms  were  generated  in  increments  of  50  nodes  up  to  250  (250  *  3  =  750  total 
nodes).  Twenty-five  instances  of  each  problem  size  were  generated  in  support  of  the 
analysis. 

The  results  in  table  3.1  include  the  problem  MPS  file  size  and  solution  time 
using  LINGO  11.0.  Twenty-five  problems  of  each  size  were  generated  and  the  mean 
and  standard  deviation  are  provided  for  the  time  required  to  generate  each  problem 
size.  The  mean  MPS  file  size  is  provided  for  each  problem  size.  Standard  deviation 
of  the  MPS  file  size  adds  no  value  to  the  analysis  and  is  omitted. 
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Nodes 

Modes 

Total  Nodes 

MPS  (MB) 

|  General 

U 

te  (sec) 
a 

Solve 

n 

(sec) 

a 

50 

3 

150 

1.3 

14.57 

0.23 

0.63 

0.18 

100 

3 

300 

5.02 

66.48 

0.82 

2.50 

0.45 

150 

3 

450 

11.24 

201.62 

2.96 

7.75 

4.58 

200 

3 

600 

19.87 

569.56 

7.87 

17.08 

8.13 

250 

3 

750 

30.94 

1364.96 

57.70 

37.48 

5.41 

Table  3.1:  Problem  Generation  and  Solution  Times 


3.5  Contributions  and  Future  Work 

This  article  has  made  two  original  contributions  to  the  state-of-the-art  in 
freight  capacity  planning.  First,  we  developed  a  new  modeling  approach  for  the 
rough-cut  capacity  planning  problem  by  modeling  this  problem  as  a  RCSP.  This 
approach  provides  a  time-scalable  method  of  capacity  planning  which  determines  a 
least  cost  allocation  of  available  capacity  to  a  network  route  which  meets  a  speci¬ 
fied  demand  (in  tons).  Secondly,  we  have  extended  this  capacity  planning  modeling 
approach  into  the  domain  of  multimodal  freight  capacity  planning.  This  extension 
allows  the  allocation  of  capacity  to  be  either  single  or  multimodal  depending  upon 
the  cost  of  such  an  allocation.  Contributions  in  this  article  demonstrate  the  validity 
of  the  proposed  modeling  method. 

Future  research  will  focus  on  generating  larger  instances  of  the  problem  to 
determine  if  exact  solutions  can  be  obtained.  If  exact  solutions  cannot  be  obtained 
in  a  reasonable  amount  of  time,  then  a  heuristic  may  need  to  be  developed  to  solve 
such  problems.  Heuristics  for  solving  the  RCSP  are  categorized  into  three  basic  types: 
path  ranking,  labeling,  or  Lagrangian  Relaxation.  Various  heuristics  use  either  pure 
forms  of  these  methods  or  a  combination  of  two  or  more  [32,42,47,59,65,82,100]. 

The  contributions  of  this  article  to  the  field  of  operations  research  could  be 
improved  by  additional  research  in  the  area  of  computer  science.  Specifically,  how  to 
quickly  generate  larger  instances  of  the  RCSP  problem  for  solution  by  commercial 
solvers.  The  limitations  experienced  during  testing  for  this  article  could  be  partially 
mitigated  by  applying  computer  science  concepts  to  aid  in  the  compact  storage  and 
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efficient  manipulation  of  the  RCSP  constraint  matrix.  However,  the  NP-complete 
classification  of  the  RCSP  will  ultimately  prevent  generating  problems  of  sufficiently 
large  size  clue  to  the  exponential  problem  growth  as  a  function  of  the  number  of 
nodes  and  arcs. 

One  area  for  future  research  is  the  extension  of  swarm  intelligence  heuristics  for 
solving  the  RCSP  and  on  the  identification  of  efficient  solutions  with  respect  to  the 
resources  of  the  RCSP.  In  the  context  of  transportation  these  resources  represent  the 
consumption  of  time  or  transportation  asset  capacity.  Identifying  efficient  solutions 
would  aid  decision  makers  in  determining  the  best  use  of  time  and  assets  to  meet 
the  goals  of  the  transportation  system. 

Additional  research  should  be  performed  in  the  areas  of  analyzing  transporta¬ 
tion  system  robustness  and  resilience  within  the  context  of  multimodal  freight  trans¬ 
portation.  Multimodal  integration  presents  a  unique  opportunity  to  reduce  redun¬ 
dancy  within  a  transportation  mode  while  maintaining  a  passive  resilience  to  meet 
demand  among  the  other  transportation  modes  within  a  system. 
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IV.  An  Ant  Colony  System  for  Solving  the  Resource  Constrained 

Shortest  Path  Problem 

4.1  Introduction 

The  Resource  Constrained  Shortest  Path  Problem  (RCSP)  is  a  variation  on 
the  classic  Shortest  Path  Problem  (SPP)  in  which  the  goal  is  to  find  a  path  of 
minimum  total  weight  connecting  a  source  node,  s,  and  a  sink  node,  t.  The  RCSP 
differs  from  the  classic  SPP  in  that  it  contains  an  additional  knapsack  constraint 
which  is  applied  to  the  sum  of  the  weights  on  the  minimum  s,t-path.  Although 
the  SPP  is  well  solved  by  algorithms  like  Dijkstra’s  algorithm,  Floyd- Warshall,  and 
the  Out-of-Kilter  algorithm,  no  efficient  algorithm  for  solving  the  RCSP  has  been 
developed  [48]. 

The  RCSP  belongs  to  the  class  of  optimization  problems  known  as  NP-Hard 
[54,110].  The  necessity  of  solving  the  RCSP  and  other,  equally  difficult  optimization 
problems  has  led  researchers  to  develop  various  metaheuristics  to  quickly  generate 
good  solutions  to  these  difficult  problems. 

The  RCSP  has  previously  been  applied  to  rough-cut  capacity  planning  in  mul¬ 
timodal  transportation  [56]  where  it  was  noted  that  sufficiently  large  instances  of 
the  RCSP  took  longer  to  generate  than  to  solve.  The  memory  required  to  generate 
the  constraint  matrix  for  the  RCSP  imposes  a  restriction  on  the  size  of  a  problem 
instance  that  can  be  generated  before  running  out  of  memory.  Even  the  largest 
instances  capable  of  being  generated  can  be  solved  to  optimality  rather  quickly. 

In  order  to  overcome  this  inherent  difficulty,  metaheuristics  can  be  used  to  solve 
the  problem.  Typically,  the  amount  of  data  required  by  a  metaheuristic  to  solve  a 
given  problem  instance  is  less  than  that  required  to  generate  the  same  problem 
constraint  matrix  for  solution  via  a  commercial  solver.  In  this  article  we  present  an 
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extension  of  the  Ant  System  (AS)  metaheuristic  that  incorporates  new  features  in 
order  to  solve  the  RCSP. 

4-1.1  Solving  the  RCSP.  Traditional  methods  for  solving  the  RCSP  have 
been  classified  into  three  categories:  Path  Ranking,  Labeling,  and  Lagrangian  Re¬ 
laxation  [47,  82] .  The  method  of  k-shortest  paths  was  developed  by  Hoffman  and 
Pavley  and  later  improved  by  Eppstein  to  run  in  0{m  +  nlogn  +  kn)  on  acyclic 
graphs  [42,59]. 

Dynamic  programming  is  one  type  of  node  labeling  method  which  operates 
recursively.  Dijkstra’s  Algorithm  is  a  classic  and  well  known  example  of  this  type 
of  solution  approach.  Irnich  and  Villeneuve  provide  a  labeling  method  for  k-cycle 
elimination  for  k  >  3  [65]. 

Lagrangian  Relaxation-based  methods  usually  begin  by  relaxing  the  resource 
constraints  which  make  the  RCSP  NP-Hard.  This  relaxation  is  solved  as  a  shortest 
path  problem  and  then  efforts  are  made  to  reduce  the  duality  gap  [32,100]. 

Many  heuristics  to  solve  the  RCSP  have  been  developed  based  on  one  or  more 
of  these  traditional  methods.  Other  solution  approaches  involve  the  use  of  meta¬ 
heuristics  like  Greedy  Randomized  Adaptive  Search  Procedure  (GRASP)  or  local 
search-based  metaheuristics  like  Evolutionary  Algorithms,  Tabu  Search,  Simulated 
Annealing,  Guided  Local  Search,  and  Iterated  Local  Search  being  extended  to  solve 
the  RCSP.  Constructive  metaheuristics  construct  feasible  solutions  and  stop  while 
local  search  metaheuristics  iteratively  improve  upon  incumbent  feasible  solutions. 

4-1.2  Constructive  and  Local  Search  Metaheuristics.  Evolutionary  Algo¬ 
rithms,  Tabu  Search,  Simulated  Annealing,  Guided  Local  Search,  and  Iterated  Local 
Search  are  examples  of  metaheuristcs  which  are  classified  as  either  constructive  or  lo¬ 
cal  search  based.  All  of  these  metaheuristics  except  GRASP,  which  is  a  constructive 
metaheuristic,  are  local  search  based. 
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4-1-3  Swarm  Intelligence  and  Ant  Colony  Optimization.  The  literature  is 
rich  with  examples  and  applications  of  biologically  inspired  metaheuristics.  Parti¬ 
cle  Swarm  Optimization,  Ant  Colony  Optimization,  and  Artificial  Bee  Colony  Algo¬ 
rithms  are  newer  examples  of  some  biologically  inspired  optimization  techniques  [22] . 

Ant  Colony  Optimization  is  a  popular  and  successful  constructive  metaheuris¬ 
tic  which  is  expressly  designed  to  solve  network  type  problems.  First  developed  in  the 
1992  dissertation  of  Marco  Dorigo,  ACO  has  been  successfully  applied  to  many  prob¬ 
lems  capable  of  being  formulated  as  shortest-path  type  problems  [36] .  The  basic  idea 
for  ACO  was  inspired  by  observing  the  behavior  of  ants  as  they  forage  out  of  the  nest 
for  a  food  source.  Each  ant,  as  it  travels  leaves  a  pheromone  trail  which  is  detected 
by  subsequent  ants.  Initially,  ants  randomly  choose  a  direction  in  which  to  move, 
then  subsequent  ants  are  influenced  by  the  pheromone  trails  left.  The  pheromones 
evaporate  over  time  and  therefore  the  pheromone  trail  on  shorter  paths  connecting 
the  nest  and  the  food  source  tends  to  be  stronger.  Eventually,  all  ants  will  choose  to 
follow  the  path  containing  the  strongest  pheromone  scent.  Empirically,  it  has  been 
observed  that  in  a  majority  of  experiments  the  shortest  path  connecting  the  nest  to 
the  food  source  is  selected  by  all  ants  after  sufficient  time  has  elapsed  [50].  ACO 
is  related  to  other  reinforcement  learning  approaches  through  artificial  pheromones 
and  evaporation  mimicking  a  process  called  stigmergy  which  makes  pheromone  trail 
strength  available  to  all  ants.  It  is  assumed  that  the  reader  is  familiar  with  the  basic 
terminology,  components,  and  operation  of  ACO  algorithms.  For  readers  requiring 
further  background  information,  the  2004  book  on  ACO  is  recommended  [37]. 

Ant  Colony  Optimization  can  include  both  constructive  and  local  search  ele¬ 
ments  and  has  been  successfully  extended  to  solve  many  shortest-path-type  prob¬ 
lems  in  the  past.  Problems  in  the  areas  of  routing,  assignment,  scheduling,  subset, 
and  machine  learning  problems  have  all  been  addressed  using  some  variation  of 
ACO  [14,17,46,68,99,114],  Many  other  articles  have  been  published  using  some 
variation  of  ACO  to  solve  a  variety  of  problems.  An  ACO  survey  paper  by  Dorigo 
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et  al.  contains  over  one  hundred  references  that  apply  AGO  variations  to  at  least 
eighteen  different  problem  types  in  five  different  problem  categories  [34], 

4-1-4  The  Basic  Ant  Colony  Optimization  Metaheuristic.  A  Metaheuristic 
is  a  heuristic  that  is  used  to  guide  other  heuristics  in  searching  a  solution  space  for  an 
optimal  solution.  A  heuristic  is  likely  to  become  “trapped”  in  local  optima,  ffence, 
metaheuristics  employ  a  variety  of  techniques  to  expand  the  search  of  the  solution 
space  and  prevent  premature  termination  of  the  search  at  local  optima.  However, 
metaheuristes  typically  do  not  guarantee  convergence  to  globally  optimal  solutions 
except  under  conditions  yielding  full  enumeration  of  the  solution  space. 

Basic  ACO  operates  by  using  a  group  of  computer  agents  (ants)  transiting  the 
solution  space  by  moving  from  a  node  to  adjacent  nodes.  The  decision  of  which 
adjacent  node  to  visit  is  determined  through  biased  random  selection.  Ants  choose 
randomly  among  remaining  nodes  but  as  the  search  progresses,  the  node  transition 
probabilities  are  biased  through  the  use  of  the  artificial  “pheromones”  deposited  by 
ants  previously  transiting  the  arc  to  a  node.  The  pheromone  deposits  are  assigned 
a  global  evaporation  rate  so  that  the  amount  of  pheromone  on  an  arc  decreases 
with  each  passing  iteration.  Obviously,  arcs  which  are  visited  more  frequently,  will 
have  higher  pheromone  concentration  as  the  search  progresses.  Ants  select  arcs  with 
greater  pheromone  concentration  with  higher  probability  than  those  arcs  with  lower 
pheromone  concentration  [37].  This  provides  the  reinforcement  aspect  of  the  search. 

Each  ant  in  the  colony  can  leverage  globally  available  information  regarding 
pheromone  concentrations.  This  indirect  method  of  communication  between  ants  in 
to  colony  is  called  “stigmergy”  and  is  the  primary  mechanism  of  the  ACO.  The  effects 
of  sitgmergy  are  studied  extensively  in  a  2005  report  prepared  for  the  Canadian 
Department  of  National  Defense  [111]. 

Several  variations  and  improvements  on  the  original  ACO  algorithm  (referred 
to  as  simple  ACO  or  S-ACO)  have  been  developed.  These  include  Ant  System  (AS), 
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Elitist  AS,  Ant-Q,  Ant  Colony  System  (ACS),  Rank-based  AS,  ANTS,  Hyper-cube 
AS,  and  Min-Max  AS  (MMAS).  Ant  System  was  initially  proposed  by  Dorigo  in 
1992  [33].  Ant  System  has  inspired  similar  metaheuristics  including  Ant-Q  and  ACS 
among  others.  Ant-Q  was  developed  by  Gambardella  and  Dorigo  in  1995  [44],  ACS 
was  inspired  by  Ant-Q  and  was  developed  in  1996  [35,45].  It  differs  from  Ant-Q  only 
by  the  initial  pheromone  value  applied  to  each  arc.  Empirical  comparisons  of  these 
ant  system  algorithms  indicate  that  the  consistently  best  performing  variants  of  the 
ant  system  algorithm  are  ACS  and  MMAS  [89]. 

ACS  is  the  most  aggressive  of  the  ACO  variations  and,  according  to  Dorigo 
generally  produces  the  best  quality  solutions  for  short  computation  times  [37].  Since 
the  method  developed  in  this  article  is  used  as  a  solution  generator  in  the  context 
of  a  larger  search  scheme,  we  have  selected  ACS  due  to  the  empirical  evidence  that 
it  produces  high  quality  solutions  in  relatively  short  periods  of  time.  Implementa¬ 
tion  details  are  discussed  in  the  next  section.  We  consider  the  RCSP  problem  as 
formulated  below: 
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N  =  {1,  2,  ...,n}  a  set  of  uniquely  labeled  nodes 

A  =  {( i,j )  :  i,  j  e  N}  a  set  of  directed  arcs  that  define 

adjacencies  for  the  nodes  in  N 
G  =  (N,  A)  is  a  graph  defined  by  the  nodes  in 

N  and  the  arcs  in  A 

R  is  the  quantity  of  resource  available 

.s  is  the  label  assigned  to  the  source 

node  where  supply  is  located 

t  is  the  label  assigned  to  the  destination 

node  where  demand  is  located 
Cij  is  the  cost  of  traversing  arc  (i,j) 

Xij  =  1  if  arc  (i,j)  is  included  in  the  path 

and  0  otherwise 


(4.1) 


Minimize  E  CijXij  (4.2a) 

s.t.  xv~  Xji  =  b(i)  VieN  (4.2b) 

j:(i,j)eA 

y.  GjXij  <  R  (4.2c) 

(i,j)eA 

b(s )  =  1,  b(t)  =  —1,  b(i )  =  0  V  i  ^  s,t  (4. 2d) 

In  this  formulation  (4.2a)  is  the  objective  function,  (4.2b)  is  the  set  of  flow 
balance  constraints,  (4.2c)  is  the  set  of  knapsack  constraints  and  (4. 2d)  are  the  flow 
forcing  constraints. 
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Merkle  et  al.  developed  an  AGO  metaheuristic  to  solve  the  Resource  Con¬ 
strained  Project  Scheduling  Problem  (RCPSP)  [84],  The  work  is  a  variation  of  the 
Ant  System  for  the  Traveling  Salesman  Problem  (AS-TSP)  originally  discussed  by 
Dorigo  and  later  by  Dorigo  et  ah  [33,36].  Improvements  to  the  AS-TSP  made  by 
Merkle  et  al.  include  combining  two  methods  of  pheromone  evaluation,  dynamic 
changes  to  the  influence  of  the  heuristic  information  on  the  probability  of  selection 
of  arcs  in  tour  construction,  and  an  option  to  “forget”  the  best  solution  found  so  far 
in  order  to  increase  exploration  during  the  search. 

Hu  et.  al.  (2010)  investigated  how  AGO  could  solve  navigation  problems  which 
account  for  multiple  psychological  expectations  of  the  driver  [63].  This  problem  is 
formulated  as  a  RCSP  and  solved  using  an  Ant  System  metaheuristic.  They  provide 
improvements  to  the  original  AS  algorithm  through  the  use  of  an  improved  transition 
probability  rule  for  the  ants  and  also  an  improved  pheromone  update  scheme  that 
consider  the  various  driver  expectations  in  route  selection. 

The  remainder  of  this  article  is  devoted  to  the  development  of  a  new  Ant 
Colony  System  heuristic  to  solve  the  RCSP  and  to  empirical  testing  of  the  heuristic’s 
performance.  An  application  to  multimodal  freight  transportation  capacity  planning 
is  presented  and  testing  is  executed  on  problem  instances  of  various  size  in  order  to 
examine  the  relationship  between  parameter  settings,  problem  size,  solution  time, 
iteration  count,  and  solution  quality.  The  article  concludes  with  a  discussion  of  the 
contributions  made  and  areas  of  future  research. 

4-2  Modeling  Approach 

4-2.1  Selecting  an  Ant  Colony  Metaheuristic.  Empirical  testing  suggests 
the  best  performing  of  the  AGO  variants  on  the  Traveling  Salesman  Problem  are  ACS 
and  MMAS  [37].  Selecting  one  of  these  variants  as  a  starting  point  for  extension  to 
solve  the  RCSP  is  based  on  comparing  the  observed  behavior  of  the  two  variants. 
Dorigo  and  Stutzle  compared  the  performance  of  ACS  and  MMAS  on  instances  of 


the  TSP  contained  in  the  TSPLIB  [37].  While  MMAS  generally  found  slightly  better 
quality  solutions  in  long  runs  of  the  algorithm,  ACS  is  the  more  aggressive  search 
strategy  and  finds  significantly  better  quality  solutions  for  short  computation  times. 
ACS  is  extended  and  applied  to  the  R.CSP. 

4  -2. 1.1  Solution  Quality  vs.  Speed.  Metaheuristics  must  be  tailored 
to  solve  specific  problems.  In  the  case  of  ACS,  there  are  parameters  that  can  be 
tuned  to  control  the  behavior  of  the  search.  Such  choices  generally  require  tradeoff 
between  computational  performance  and  solution  quality.  Some  major  configuration 
decisions  in  implementation  are  whether  or  not  to  implement  a  local  search  strategy 
and  if  local  search  is  used  then  determining  pheromone  update  type  (Darwinian  vs. 
Lamarckian),  use  of  data  structures,  and  type  of  heuristic  information  about  network 
arcs. 

Implementation  of  a  local  search  begins  with  a  tour  constructed  by  an  ant 
(si)  and  through  the  local  search  process  yields  an  improved  path  (52).  Reinforcing 
the  pheromone  trail  on  .sq  is  referred  to  as  a  Darwinian  pheromone  update  while 
reinforcing  the  pheromone  trail  along  the  arcs  corresponding  to  S2  is  a  Lamarckian 
update. 

In  this  research  no  local  search  is  implemented  following  initial  path  construc¬ 
tion.  Eliminating  local  search  produces  a  time  savings  of  about  ten  percent  [37]. 
Heuristic  information  is  used  in  tour  construction.  The  method  of  tour  construction 
used  in  ACS  determines  for  each  ant  k  at  city  i  the  city  j  that  is  visited  next  by: 


I  argmaxleNk  {[Tii]a[r}ii]0}  ,  if  q  <  <?o;  ,,  oN 

J  =  <  *  (4-3) 

I  J,  otherwise 

In  4.3,  the  parameter  r/h  is  the  multiplicative  inverse  of  the  Euclidean  dis¬ 
tance  between  i  and  j,  (r/h  =  1  /dij).  This  Euclidean  distance  heuristic  provides  a 
more  aggressive  search  strategy  than  simply  using  the  multiplicative  inverse  of  arc 
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distance  to  provide  the  heuristic  information.  The  parameter  t-  is  the  strength  of 
the  pheromone  on  arc  N A  is  the  neighborhood  of  ant  k  while  at  node  i,  q  is 

a  uniformly  distributed  random  variable  on  the  interval  [0,1],  (0  <  go  <  1),  f3  is 
the  tuning  parameter  for  the  weight  of  the  heuristic  value  relative  to  the  pheromone 
strength,  and  a  is  the  pheromone  weight  parameter.  Also  let  J  be  a  random  variable 
from  the  probability  distribution: 


Y2l£Nk  ['TilY'VlilY 


,  if  jeN* 


(4.4) 


Pheromone  levels  are  updated  after  each  ant  has  constructed  a  path.  The  two 
update  operations  applied  to  the  pheromones  are  evaporation  and  deposit.  Collec¬ 
tively,  these  two  operations  provide  the  stigmergistic  aspects  of  the  heuristic.  Stig- 
mergy  is  manipulated  by  changing  the  rate  and  method  for  evaporation  and  deposit. 
All  ants  have  access  to  this  common  information  and  therefore  communicate  indi¬ 
rectly  through  updates  to  these  stigmergic  parameters.  Evaporation  controls  how 
quickly  the  corporate  “memory”  of  estimated  arc  quality  is  “forgotten.”  Deposit  up¬ 
dates,  discussed  in  the  following  section,  provide  a  means  of  updating  the  collective 
corporate  memory. 


Finally,  ASC  differs  form  other  AGO  metaheuristics  in  three  basic  ways.  First, 
it  uses  a  more  aggressive  search  strategy  in  the  construction  step  by  heavily  exploiting 
the  learned  knowledge  of  other  ants  through  the  pheromone  trails.  Second,  ACS  only 
employs  pheromone  deposit /evaporation  on  the  best  path  found  so  far  in  the  search. 
Finally,  pheromone  strength  is  reduced  on  an  arc  immediately  after  it  is  traversed 
by  an  ant  which  drives  exploration  by  making  previously  explored  arcs  less  likely  to 
be  selected. 


In  this  research,  we  use  a  designed  experiment  to  optimize  the  parameter  set¬ 
tings  used  in  ACS  as  applied  to  the  RCSP.  Typical  implementations  of  ACS  assume 
a  value  of  a  =  1.  In  our  experimentation  we  treat  a  as  a  variable  value  to  which 
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optimization  is  applied  along  with  the  other  search  parameters  in  order  to  determine 
the  best  setting  of  all  parameters  relative  to  the  RCSP. 

4-2.2  Extension  of  ACS  to  Solve  RCSP  Problem.  Ant  Colony  Metaheuris¬ 
tics  are  often  used  to  construct  Hamiltonian  Circuits  in  solving  the  traveling  salesman 
problem.  Such  a  search  explicitly  requires  that  each  node  in  the  network  is  visited 
exactly  one  time.  In  constructing  shortest  paths  it  is  only  necessary  to  specify  a 
source  and  sink  node  and  perform  the  search  until  the  sink  node  has  been  visited 
at  which  point  the  search  terminates  and  the  current  path  is  returned  as  a  feasible 
solution. 

In  ACS,  both  a  global  and  a  local  pheromone  update  are  used.  The  local 
pheromone  update  is  done  by  all  ants  immediately  after  traversing  an  arc.  The 
global  pheromone  update  is  done  once  per  iteration  only  by  the  ant  that  has  found 
the  best  path  so  far. 

An  iteration  consists  of  all  operations  required  for  each  ant  in  the  colony  to 
build  a  path,  determine  the  best  path  from  among  these  paths  ( sa ),  compare  (and 
replace  if  necessary)  to  the  best  path  so  far  over  all  iterations  (sbs),  and  complete 
the  local  and  global  pheromone  updates,  both  deposits  and  evaporation  procedures. 

Global  pheromone  updates  are  applied  at  the  end  of  each  iteration  according 
to  the  following  replacement  operation:  <—  (i  —  p)Tij  +  pArff ,  V(i,  j)  E  Pbs.  In  this 

operation,  A rbf  =  l/(Cbs),  where  Cbs  is  the  length  of  the  best  path  found  so  far,  Pbs 
is  the  best  path  so  far,  and  p  E  (0, 1]  is  the  pheromone  evaporation  rate  parameter. 
Local  pheromone  updates  are  applied  according  to:  rtJ  E-  (1  —  f)rij  +  £ro.  In  this 
operation  £  E  (0,1]  is  an  evaporation  parameter.  The  local  pheromone  update  is 
applied  to  an  arc  immediately  after  an  ant  crosses  the  arc.  The  effect  of  the  local 
pheromone  update  is  to  decrease  the  desirability  of  an  arc  after  it  is  traversed  in 
order  to  drive  diversification  of  the  search  and  avoid  search  stagnation.  The  effect  of 
the  global  pheromone  update  is  to  reinforce  the  pheromone  on  arcs  that  are  on  the 
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best  path  found  so  far  in  the  search.  The  local  pheromone  update  prevents  search 
stagnation  by  encouraging  exploration  of  previously  unvisited  arcs. 

Since  ACS  is  usually  applied  to  the  TSP,  the  r0  parameter  is  typically  set  to 
a  value  of  1  /(riCnn)  where  n  is  the  number  of  cities  in  the  TSP  and  Cnn  is  the 
length  of  a  tour  obtained  using  a  “nearest  neighbor”  heuristic  [37].  In  shortest-path 
problems  we  do  not  have  an  a  priori  value  for  either  n  or  Cnn  and  will  have  to 
estimate  these  parameters.  Options  for  estimating  n  might  include  shortest  path 
heuristics  like  A*  or  solving  the  shortest  path  problem  directly  using  a  shortest 
path  algorithm  like  Dijkstra’s  Algorithm.  Either  of  these  options  could  significantly 
increase  computational  time. 

Estimates  for  n  and  Cnn  are  obtained  by  approximating  the  number  of  “hops” 
on  the  shortest  path  (nst),  and  the  length  of  the  shortest  s,t-path  estimated  as  the 
Euclidean  s,t  distance  ( Est ).  Now  our  estimate  is  taken  as:  r0  =  l/(nstEst).  nst  is 
obtained  by  dividing  Est  by  the  average  length  of  an  arc  in  the  network  (a/).  To 
obtain  a;  we  divide  the  sum  of  all  network  arc  distances  (t)  by  the  number  of  arcs 
in  the  network  (m).  We  can  calculate  To  directly  for  any  network  encountered  as: 
r0  =  t  /  (m(Est)2) . 

Dorigo  and  Stutzlc  suggest  using  the  multiplicative  inverse  of  arc  distance, 
rjij  =  j-,  as  a  heuristic  measure  of  the  desirability  of  arc  (i,j)  [37].  If  we  assume 
Euclidean  distance  then:  disij  =  [(x*  —  Xj)2  +  (yt  —  yj )2]1/2.  Another  heuristic  is: 


Va 


i 

disj  d 

1 


if  disjd  ^  0 
if  disjd  =  0 


(4.5) 


Where  disjd  is  the  Euclidean  distance  from  node  j  to  the  destination  node  d. 
This  Euclidean  distance  heuristic  for  use  in  AGO  was  first  proposed  by  Hu  et.  al.  [63] 
to  make  the  search  less  “myopic”  at  each  iteration  by  encouraging  selection  of  nodes 
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that  are  closer  to  the  destination  node.  The  relative  importance  of  the  pheromone 
and  heuristic  values  are  mediated  by  the  parameters  a  and  /3. 

Traditional  Ant  Colony  Optimization  constructs  feasible  paths  by  accounting 
for  node  adjacency  in  non-complete  graphs.  The  RCSP  has  the  knapsack  constraint 
that  imposes  an  additional  restriction  on  path  feasibility.  A  given  path  is  feasible 
only  if  it  is  also  resource-feasible.  Therefore,  the  neighborhood  of  the  current  node 
consists  of  nodes  that  are  adjacent  and  resource-feasible. 

Two  mechanisms  are  used  to  ensure  that  resource-feasible  paths  are  guaranteed 
to  be  constructed  if  they  exist.  First,  resource  consumption  on  partially  constructed 
paths  is  maintained  for  each  ant,  making  possible  the  identification  of  non- resource- 
feasible  nodes  in  each  neighborhood.  Secondly,  backtracking  is  implemented  to  allow 
ants  to  move  to  a  previously  visited  node  if  no  resource-feasible  nodes  exist  in  the 
current  neighborhood.  This  approach  guarantees  construction  of  resource  feasible 
paths  if  such  paths  exist. 

4-3  Testing  the  ACS  Metaheuristic 

Our  ACS  metaheuristic  was  tested  in  three  stages.  The  first  stage  involves  de¬ 
termining  the  parameter  values  for  the  ACS.  The  second  stage  involves  running  ACS 
on  problems  with  known  optimal  solutions  to  characterize  the  quality  of  the  ACS 
solutions  along  with  the  speed  with  which  they  are  obtained.  Finally,  ACS  is  tested 
on  large  problem  instances  to  demonstrate  that  the  ACS  can  quickly  find  feasible 
solutions  to  problems  that  cannot  be  generated  for  solution  using  a  deterministic 
binary  integer  programming  approach. 

4-3.1  Selection  of  Parameter  Values.  An  ACS,  offers  two  primary  objec¬ 
tives:  “ running  time”  and  “ percent  above  the  minimum  objective  function  value.” 
These  response  variables  are  affected  by  the  problem  and  by  the  ACS  parameter 
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Parameter  Name 

Description 

Limits  (Lower, Upper)* 

a 

controls  the  relative  influ¬ 
ence  of  the  pheromone  trail 

[0,1,1] 

P 

controls  the  relative  influ¬ 
ence  of  the  heuristic  infor¬ 
mation 

[1,6] 

max  iterations 

the  number  of  complete 
search  iterations  each  ant 
colony  should  perform 

[10, 1000] 

num  ants 

the  size  of  the  ant  colony 
used  in  the  search 

[5, 1000] 

num  neighbors 

the  number  of  nearest 
neighbors  to  consider  in 
selecting  each  successive 
node 

[5,25] 

% 

pseudorandom  proportional 
action  choice  rule 

[0,1,1] 

P 

pheromone  evaporation  rate 

[0.1,  0.9] 

z 

local  pheromone  trail  up¬ 
date  rule 

[0,1,1] 

*  limits  as  suggested  by  Dorigo  and  Stutzle  [37] 


Table  4.1:  Parameter  Value  Ranges 
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values.  Design  of  Experiments  (DOE)  and  Response  Surface  Methodology  (RSM) 
provide  tools  to  find  input  parameter  values  [92], 

Optimizing  multiple  response  variables  can  be  accomplished  by  fitting  models 
to  each  response  variable  and  then  jointly  optimizing  both  responses  using  desirabil¬ 
ity  functions  [31].  Since  ACS  has  stochastic  components,  replications  are  used.  The 
focus  of  of  our  initial  experimentation  is  to  characterize  the  solution  space.  Rather 
than  using  a  28  full  factorial  design  of  256  experimental  runs,  we  opt  for  a  central 
composite  design  (CCD)  [92],  This  design  ensures  independent  estimates  of  all  main 
effects  and  two-way  interaction  terms.  However,  some  of  our  experimental  parame¬ 
ters  must  be  set  to  integer  values  and  other  parameters  are  being  experimented  with 
at  the  their  operating  extremes.  Since  the  CCD  axial  points  may  not  correspond  to 
integer  parameter  values  we  opt  to  use  a  non-rotatable  variation  of  the  CCD  called 
a  face-centered  cube  design  (FCD).  Although  the  FCD  is  non-rotatable,  it  has  fairly 
stable  prediction  variance  throughout  the  design  region  when  two  center  runs  are 
used.  No  significant  improvement  in  prediction  variance  is  gained  beyond  two  center 
runs  and  therefore  we  elect  to  use  exactly  two  center  runs  in  the  FCD.  We  have 
selected  a  fractional  FCD  composed  of  a  28_1  with  128  runs,  18  axial  points,  and 
two  center  points.  The  total  number  of  experiments  is  146. 

Empirical  testing  by  Hartlage  and  Weir  indicates  that  the  largest  problem 
instance  able  to  be  generated  for  solution  by  a  commercial  solver  had  750  nodes  [56]. 
All  attempts  to  generate  larger  problems  failed  due  to  memory  issues.  We  ran  the 
ACS  experiments  here  using  a  randomly  generated  network  with  750  total  nodes  and 
ninety  percent  arc  density.  The  table  below  shows  the  fit  statistics  and  the  analysis 
of  variance  information  obtained  using  JMP. 


Not  surprisingly,  the  model  fit  for  for  the  response  of  “percent  over  optimal”  is 
weak.  This  is  driven  mainly  by  the  propensity  of  the  ACS  metaheuristic  to  converge 
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R2 

0.375135 

JD  2 

adjusted 

0.102916 

RMSE 

7.070444 

Observations 

146 

Table  4.2:  Fit  Statisitcs  for  %  over  optimal 


Source 

DF 

SS 

MS 

F  ratio 

Model 

44 

3031.2051 

68.8910 

1.3781 

Error 

101 

5049.1084 

49.9912 

Prob  >  F 

Total 

145 

8080.3135 

0.0953 

Table  4.3:  ANOVA  Table  for  %  over  optimal 


to  an  optimal  solution  in  an  overwhelming  number  of  test  instances,  even  when 
varying  the  parameter  settings.  This  tendency  indicates  that  the  ACS  metaheuristic 
is  fairly  robust  to  parameter  values.  The  plot  of  predicted  versus  actual  values 
provides  visual  evidence  of  this  assertion.  Notice  the  appearance  of  a  horizontal 
grouping  of  test  points  whose  actual  value  was  “zero.”  That  is,  these  points  are 
test  points  for  which  the  predicted  value  was  greater  than  the  actual  value.  Stated 
another  way,  the  metaheuristic  outperformed  the  predicted  performance  in  these 
cases.  Figure  4.1  provides  a  concise  graphical  depiction  of  this  relationship. 


Figure  4.1:  Percent  over  optimal  for  actual  vs.  predicted 
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The  table  below  summarizes  this  empirical  result.  Notice  that  over  sixty- 
percent  of  the  time,  regardless  of  the  parameter  combination  used  in  the  experiment, 
ACS  found  a  solution  that  is  within  just  four  percent  of  the  optimal  solution. 


optimal  runs 

49  of  146 

<  0.01 

86  of  146 

<  0.04 

88  of  146 

Table  4.4:  ACS  Performance  during  Initial  Experimentation 

Although  the  relationship  between  predicted  and  actual  values  is  statistically 
weak,  the  p- value  indicates  that  with  an  a  of  0.1  the  model  still  detects  variability 
caused  by  the  different  parameter  settings. 

The  fit  for  the  response  of  “run  time”  was  significantly  better  than  the  model 
for  “percent  over  optimal.”  Obviously,  the  run  time  of  the  heuristic  is  primarily  a 
function  of  many  controllable  factors  within  ACS.  The  size  of  the  ant  colony,  the 
number  of  iterations  performed,  and  the  number  of  neighbors  in  the  restricted  nearest 
neighbor  list  for  each  node  are  obvious  drivers  of  this  response.  All  values  for  run 
times  are  in  units  of  seconds. 


The  high  adjusted  R2  value  and  the  small  p- value  indicate  that  the  changes 
in  parameter  values  significantly  affect  the  variability  in  running  time,  and  that  the 
fitted  model  shows  the  strong  relationship  between  predicted  and  actual  values. 

The  parameter  settings  that  yield  maximum  desirability  of  the  two  responses, 
are  presented  in  the  following  table.  These  values  will  be  used  in  all  remaining 
experimentation. 


R2 

0.911130 

TD  2 

adjusted 

0.872414 

RMSE 

42.0135 

Observations 

146 

Table  4.5:  Fit  Statistics  for  Run  Time 
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Source 

DF 

SS 

MS 

F  ratio 

Model 

44 

1827781.7 

41540.5 

23.5339 

Error 

101 

178278.5 

1765.1 

Prob  >  F 

Total 

145 

2006060.2 

<0.0001 

Table  4.6:  ANOVA  Table  for  Run  Time 


Figure  4.2:  Run  Time  for  actual  vs.  predicted 


Coded  Value 

Natural  Value 

a 

0.1873619 

0.6343 

0 

-1 

1 

Maxlterations 

-0.698048 

7  (rounded) 

NumAnts 

-1 

5 

NumNeighbors 

-1 

2 

<?o 

0.4758798 

0.7641 

P 

1 

0.9 

£ 

1 

1 

Table  4.7:  Final  ACS  Parameters  Maximizing  Desirability 


The  next  section  presents  results  obtained  from  testing  the  ACS  metaheuristic 
using  the  parameter  settings  in  the  table  above  which  are  found  to  jointly  maximize 
the  desirability  of  the  ACS  metaheuristic  running  time  and  solution  quality,  the  two 
responses  of  interest. 


4-4  Results 

The  results  in  this  section  are  based  on  two  types  of  testing.  The  goal  of  the  first 
group  of  tests  is  to  provide  empirical  verification  of  the  final  ACS  parameter  settings 
affect  on  the  two  responses.  Certain  combinations  of  parameter  settings  produce 
significantly  longer  running  times  even  on  smaller  problems  able  to  be  solved  to 
optimality  by  solver  software. 

Due  to  the  size  of  the  constraint  matrix  generated  for  relatively  small  problem 
instances,  not  only  was  the  process  of  problem  generation  time  consuming,  but  no 
constraint  matrix  for  BIP  exceeding  250  *  3  =  750  total  nodes  was  successfully 
generated  due  to  memory  overrun  issues.  Due  to  this  limitation,  testing  ACS  solution 
quality  required  maintaining  small  problem  sizes  (750  total  nodes)  in  order  to  obtain 
optimal  solutions  for  comparison. 

4-4-1  Small  Example  Problems.  Using  the  settings  obtained  through  the 
method  of  maximum  desirability  we’ll  now  run  ACS  against  LINGO  to  compare 
speed  and  solution  quality  for  different  problems.  Density  of  arcs  refers  to  the  number 
of  arcs  in  the  network  relative  to  the  number  of  arcs  contained  in  a  complete  graph 
on  the  same  number  of  nodes.  We  use  an  arc  density  of  ninety  percent  for  purposes  of 
experimentation.  This  density  was  selected  through  trial  and  error.  Greater  densities 
tend  to  produce  single  arc  solutions  which  are  rather  uninteresting  as  they  relate 
to  multimodal  paths  while  less  dense  solutions  sometimes  result  in  disconnected 
networks  with  no  feasible  solutions.  In  real-world  problems,  connectivity  would  be 
less  of  an  issue  since  the  input  to  the  problem  is  typically  a  transportation  network 
currently  in  use.  For  a  complete  discussion  of  graph  completeness  and  density  see 
the  text  by  West  [12].  Table  4.8  summarizes  the  results  of  this  testing  by  providing 
a  summary  of  the  running  time  (in  seconds)  for  ACS  and  BIP  along  with  a  column 
showing  the  performance  of  ACS  in  terms  of  solution  quality  (given  in  terms  of 
percent  over  optimal). 
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For  this  experiment  test  we  perform  twenty-five  replications  using  randomly 
generated  networks  containing  250  nodes  and  three  modes.  All  networks  have  90 
percent  arc  density.  The  results  in  Table  4.8  indicate  that  the  parameter  values 
selected  by  the  method  of  maximum  desirability  produce  solutions  are  superior  to 
those  obtained  in  the  parameter  selection  experiments  in  both  running  time  and 
solution  quality. 


Parameter  Select  Expt 

Max  Des  Settings  Expt 

Num  runs 

146 

25 

Number  of  Optimal  Runs 

49  =  33.56% 

12  =  48% 

Runs  within  1%  of  Opt 

86  =  58.90% 

14  =  56% 

Runs  within  3%  of  Opt 

88  =  60.27% 

20  =  80% 

Worst  run  Pers  over  opt 

3674.71% 

11% 

Mean  Run  Time  (sec) 

65.49 

0.9974 

StDev  Run  time 

117.62 

0.5732 

99%  Cl  for  Run  Time 

(40.38,90.60) 

(0.70,1.29) 

Mean  Percent  over  Opt 

338.97% 

2.00% 

StDev  Percent  over  Opt 

746.50% 

3.16% 

99%  Cl  for  Percent  over  Opt 

(179.58,498.36) 

(0.37,3.63) 

Table  4.8:  Max  Desirability  vs.  Parameter  Select  Expt  Settings 


Note  that  the  99%  confidence  intervals  on  the  mean  run  time  and  percent  over 
optimal  for  the  two  experiments  do  not  overlap  indicating  statistical  superiority  of 
the  max  desirability  parameter  settings. 

The  limiting  factor  in  the  use  of  deterministic  approaches  to  the  RCSP  is  in 
generating  the  constraint  matrix  for  the  problem.  The  previous  tests  demonstrate 
empirically  that  ACS  produces  solutions  (even  removing  generation  time  from  con¬ 
sideration)  more  quickly  than  deterministic  approaches  for  small  problems.  The 
remaining  set  of  experiments  is  intended  to  demonstrate  that  ACS  is  able  to  quickly 
generate  and  solve  problems  that  are  too  large  to  be  generated  for  solution  via  BIP. 

4-4-%  Larger  Example  Problems.  We  begin  this  section  by  solving  problems 
with  three  transportation  modes  and  600  nodes  for  a  total  f  1800  nodes.  Nodes  are 
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progressively  incremented  by  25  *  3  =  75  nodes  for  each  successive  experiment. 
Twenty- five  problems  of  each  size  with  90%  arc  density  are  solved  and  the  results 
are  presented  in  Table  4.9.  These  larger  probelms  demonstrate  the  applicability  of 
this  solution  approach  for  realistically-sized  problem  instances. 


Total  Nodes 

Gen  (sec) 

Gen  <x 

Solve  (sec) 

Solve  a 

1800 

2.2290 

0.0822 

3.2798 

0.9335 

1875 

2.3810 

0.0750 

3.3292 

0.4862 

1950 

2.5576 

0.0705 

4.4426 

0.6586 

2025 

2.7453 

0.1385 

3.6931 

0.4839 

2100 

3.0298 

0.0573 

4.5830 

0.4793 

2175 

3.1771 

0.1266 

4.5443 

0.5391 

2250 

3.4552 

0.1202 

4.5717 

0.5026 

2325 

3.6262 

0.0870 

6.0689 

0.7672 

2400 

3.8958 

0.1316 

5.7184 

0.5939 

2475 

4.1083 

0.1147 

5.1899 

0.4678 

2550 

4.4337 

0.1224 

5.1891 

0.5870 

2625 

4.6430 

0.1288 

6.1895 

0.4923 

2700 

4.9430 

0.1186 

8.6282 

0.6248 

2775 

5.2524 

0.1430 

7.6742 

0.6755 

2850 

5.5732 

0.1367 

7.0104 

0.4993 

2925 

6.0586 

0.1643 

7.8940 

0.6168 

3000 

6.1649 

0.1722 

7.2732 

0.6885 

3075 

6.4917 

0.1337 

8.6660 

0.6033 

3150 

6.7733 

0.1803 

8.0548 

1.1692 

3225 

7.0783 

0.1586 

8.8562 

0.6442 

3300 

7.4648 

0.1516 

9.3473 

0.5789 

3375 

7.8142 

0.2569 

9.4269 

0.7746 

3450 

8.1275 

0.3497 

11.1530 

0.8501 

3525 

8.3575 

0.2204 

11.3913 

0.9032 

3600 

8.8417 

0.3062 

10.8643 

0.6984 

Table  4.9:  Mean  and  a  for  Large  Experiments 


4-5  Conclusions  and  Future  Work 

This  research  developed  an  extension  of  the  Ant  Colony  System  metaheuristic 
for  solving  the  Resource  Constrained  Shortest  Path  Problem.  Empirical  results  show 
strong  potential  for  using  this  metaheuristic  to  quickly  obtain  high-quality  solutions 
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to  realistic  sized  problem  instances  which  are  too  large  for  solution  by  deterministic 
methods. 

Future  research  should  focus  on  two  main  areas.  The  first  is  on  implementing 
the  ACS  heuristic  in  a  faster  programming  language  like  C++.  Gains  in  speed  for 
ACS  are  nearly  guaranteed  by  implementing  it  in  a  compiled  rather  than  an  inter¬ 
preted  language.  Secondly,  further  research  into  compact  problem  representation 
would  be  useful  since  it  is  problem  generation  due  to  computer  memory  overruns 
rather  than  solution  speed  that  limits  the  problem  size  able  to  be  solved. 
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V.  Finding  Near  Efficient  Solutions  to  Lexicographic, 

Multi- Objective,  Resource  Constrained  Shortest  Path  Problems 

5. 1  Introduction 

Capacity  planning  in  transportation  is  an  integral  and  important  function  for 
mid-range  planning  or  fleet  utilization.  Rough-Cut  Capacity  Planning  (RCCP)  in 
multimodal  freight  transportation  has  previously  been  addressed  by  Hartlage  and 
Weir  and  is  solved  using  either  deterministic  binary  integer  programming  solvers 
or  metaheuristics  [56].  Since  the  RCSP  is  known  to  be  NP-Hard,  sufficiently  large 
problem  instances  require  fast  running  metaheuristics  to  quickly  provide  high-quality 
solutions.  A  new  metaheuristic,  ACS-RCSP,  based  on  the  Ant  Colony  System  (ACS) 
traditionally  used  to  solve  the  Traveling  Salesman  Problem  (TSP)  was  developed  by 
Hartlage  and  Weir  [55]. 

In  RCCP  for  multimodal  freight,  three  controllable  parameters  ultimately  de¬ 
termine  resource  requirements  and  constraint  right-hand-sides  (RHS).  These  param¬ 
eters  are: 

•  number  of  tons  (in  millions)  of  freight  to  transport, 

•  number  of  transportation  assets  available  in  each  transportation  mode,  and 

•  length  (in  days)  of  the  planning  horizon. 

We  use  the  term  “controllable”  to  distinguish  between  parameters  that  directly 
affect  the  RHS  (or  arc  weights)  and  are  completely  determined  by  the  decision  maker, 
and  those  parameters  that  affect  the  RHS  but  are  primarily  determined  by  factors  not 
directly  under  the  control  of  the  decision  maker.  Parameters  in  the  latter  category 
are  referred  to  as  “uncontrollable.”  Some  examples  of  parameters  in  this  category 
include: 

•  productivity  of  each  transportation  asset  type, 
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•  utilization  rate  of  each  transportation  asset  type,  and 

•  block  speed  of  each  transportation  asset  type. 


To  illustrate  these  parameters,  consider  asset  utilization  or  “UTE  rate.”  The 
UTE  rate  for  aircraft  is  defined  by  Air  Mobility  Command’s  Air  Mobility  School  as 
“the  total  hours  of  capability  a  fleet  of  airlift  aircraft  can  produce  in  a  day  expressed 
in  terms  of  per  primary  authorized  aircraft.”  Similar  definitions  are  available  for 
other  modes  of  transportation. 

Note  that  higher  UTE  rates  are  not  necessarily  “better”  in  any  sense  of  the 
word.  A  decision  maker  may  elect  to  increase  UTE  to  an  unsustainable  “surge”  UTE 
rate  for  a  short  period  of  time  but  this  is  a  strictly  short-term  situation. 


UTE  rate  is  affected  by  many  factors.  A  list  of  twenty  four  different  factors 
affecting  aircraft  UTE  rate  illustrates  that  some  of  the  factors  influencing  UTE  rate 
are  controllable  and  some  are  not.  This  list  is  provided  for  ease  of  reference  and  in 
order  to  clearly  illustrate  the  point  that  UTE  rate  is  not  directly  controllable: 


Average  Ground  Time 
Airspeed  En  Route 
Crew  Availability 
Crew  Stage  Base  policies 
Reserve  Call-Up  Schedule 
Maintenance  Manpower 
Ramp  Space  Constraints 
JCS  Withhold  Levels 
Air  Traffic  Control  Delays 
Weather 

Mission  Capable  Rate 
Aircraft  Maintainability 
Aircraft  Generation  Schedule 


Average  Mission  Time  and  Leg  Length 
Crew  Ratio 

Crew  Augmentation  policies 
Active  and  Reserve  Force  Mix 
Spares  ande  Resupply  Availability 
Scenario  Resource  Constraints 
PAA  Airframes 
Aerial  Refueling  Policies 
Political  and  diplomatic  clearance  delays 
Airfield  operating  hours 
Aircraft  Reliability 
Aircraft  Availability 
Sch.  Maint.  Interval  and  Duration 


Table  5.1:  Table  of  LITE  Factors 


5.1.1  Notions  of  Optimality.  Multicriteria  optimization  problems  define 
optimality  in  several  ways.  We  must  first  define  how  the  various  objective  compo¬ 
nents  are  to  be  compared  before  a  notion  of  optimality  or  efficiency  can  be  applied. 
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Objectives  are  compared  explicitly  by  ordering  them  according  to  some  set  of  crite¬ 
ria.  For  a  complete  discussion  on  ordering  objectives  and  types  of  orders,  refer  to 
Ehrgott  [40]. 

We  are  concerned  with  different  lexicographic  orderings  of  the  objective  com¬ 
ponents.  Lexicographic  ordering  assigns  a  priority  to  each  objective  component  so 
that  a  small  increase  in  a  higher  priority  component  outweighs  a  large  increase  in  a 
lower  priority  component.  Examining  the  efficient  solutions  that  result  from  different 
objective  component  prioritizations  is  a  major  item  if  interest  in  this  research. 

Depending  upon  the  size  of  the  RCSP  problem  instance  (number  of  network 
nodes  and  density  of  network  arcs),  it  may  or  may  not  be  practical  to  locate  these 
efficient  solutions  using  a  binary  integer  programming  (B1P)  approach.  In  cases 
where  the  problem  instance  is  too  large  for  BIP  to  be  practical  a  metaheuristic,  such 
as  the  ACS- RCSP  may  be  used  to  generate  near-optimal  solutions  [55].  In  cases 
where  a  suboptimal  solution  is  returned  by  ACS-RCSP,  we  consider  such  solutions 
to  be  “near-efficient” . 

5.1.2  Problem  Statement.  The  problem  under  investigation  in  this  research 
is  the  location  and  identification  of  efficient  solutions  to  the  RCSP  for  multimodal 
freight  transportation.  Locating  such  solutions  requires  defining  some  notion  of 
optimality.  In  the  event  that  a  lexicographic  ordering  of  the  transportation  asset 
types  is  provided,  near-efficient  solutions  are  generated  with  respect  to  the  ordering. 
If  no  objective  ordering  is  provided,  the  near-efficient  solutions  are  generated  with 
respect  to  minimizing  costs  by  considering  the  price,  in  dollars,  per  million-ton- 
rnile  for  each  transportation  mode,  and  as  such  an  implicit  lexicographic  objective 
ordering  is  used. 

5.1.3  Review  of  the  Literature. 
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5. 1.3.1  Multimodal  Freight  Capacity  Planning.  Capacity  planning 
in  multimodal  freight  transportation  has  received  relatively  little  attention  in  the 
literature  to  date  [79].  In  fact,  much  of  the  previous  research  addressing  multiple 
transportation  modes  focuses  on  mode  selection  rather  than  integration.  A  recent 
survey  of  the  transportation  mode  choice  and  carrier  selection  literature  reviewed 
over  ninety  separate  references  spanning  a  period  of  over  forty  years  [83]! 

The  2005  article  by  Park  and  Regan  developed  a  bi-level  programming  for¬ 
mulation  for  multimodal  freight  capacity  planning.  However  the  approach  uses  an 
“average  load  factor”  for  each  transportation  mode  and  is  therefore  not  additive  with 
respect  to  capacity  over  the  entire  fleet  [96]. 

In  2006  Sun  et  al.  use  an  existing  linear  programming  model  called  MAX- 
CAP  which  seeks  to  maximize  system  capacity  for  multiple  origin-destination  pairs 
subject  to  resource  and  capacity  constraints.  They  analyze  system  “flexibility”  and 
performance  under  degradation  but  consider  only  a  single  mode  of  transportation, 
movement  by  rail  [104], 

Morlok  and  Riddle  address  the  maximization  of  multimodal  system  capacity 
by  defining  “traffic  lanes”  that  represent  aggregated  paths  and  modes.  Aggregation 
is  done  based  on  origin/destination  pairs  that  are  geographically  close  to  each  other. 
Capacity  of  the  system  is  estimated  by  summing  over  the  aggregated  traffic  lanes  [88] . 

Unnikrishnan  and  Waller  analyze  the  effect  of  of  demand  and  capacity  uncer¬ 
tainty  on  rail  network  capacity,  seeking  to  maximize  the  fraction  of  demand  satisfied 
for  each  origin/destination  pair  in  the  network  [106]. 

Cambridge  Systematics  Inc.  published  a  manual  that  discusses  multimodal  net¬ 
works  for  both  freight  and  passenger  transportation.  This  manual  uses  multimodal 
“corridors”  which  are  essentially  channel  routes  that  utilize  multiple  transportation 
modes.  Capacity  estimation  of  the  network  is  done  by  analyzing  the  corridor  capacity 
by  analysis  of  individual  link  capacity  in  each  multimodal  corridor  [105]. 


76 


Andersen  and  Crainic  address  the  design  and  optimization  of  a  multimodal 
scheduled  service  network.  The  approach  used  in  this  article  is  to  maximize  system 
throughput  by  coordinating  and  optimizing  the  departure  times  for  multiple  vehi¬ 
cle  fleets  utilizing  a  common  network  [4].  No  literature  source  reviewed  addressed 
identification  of  efficient  solutions  under  different  parameter  value  assumptions. 

5. 1.3. 2  RCSP  Classification  and  Solution.  The  basic  RCSP,  formu¬ 
lated  below,  is  NP-Complete  as  shown  by  Handler  and  Zang  [54],  To  simplify  the 
presentation  of  the  single  mode  resource  constrained  shortest  path  formulation  the 
following  parameters  and  vectors  are  defined:  N  =  {1,2 ,  ...,n}  is  the  set  of  nodes 
representing  shipping  locations,  C  =  {( l,v )  :  (l,v)  is  an  arc  in  the  graph}  is  the  set 
of  feasible  directed  arcs  defined  over  the  elements  in  N,  G  =  (N,  C )  is  the  graph 
consisting  of  nodes  in  N  and  edges  in  C ,  and  E  is  the  n  x  \C\  node-arc  adjacency 
matrix  where  each  column  of  E  has  a  1  in  the  Ith  row  and  a  —1  in  the  Vth  row  corre¬ 
sponding  to  adjacency  matrix  element  Eiv.  Now  define  Eiv  along  with  the  following 
sets,  matrices,  and  scalars  as: 
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Eij  - 


1  if  l  G  N  is  the  tail  of  arc  j  G  C 

—  1  if  Z  G  N  is  the  head  of  arc  j  e  C  V  l  G  iV,  j  G  C 

0  otherwise 


1  if  ESj  =  1 
0  otherwise 


V  iec 


g-f(j)  = 


1  if  Etj  =  -1 

0  otherwise 


V  jeC 


Q  = 


(5.1) 


h 


the  cost  (in  dollars)  per  MTM 


r  : 


the  number  of  MTM  available 


B_s  :  the  1  x  \C\  vector  of  arc  weights 


x  :  the  1  x  \C\  vector  of  decision  variables 


J  = 


B_s 

~h~ 


The  single  mode  formulation  can  now  be  expressed  as: 
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Minimized  s  x 


(5.2a) 

subject  to :Ex  =  0  (5.2b) 

Qx  =  1  (5.2c) 

Jx  <  r  (5. 2d) 

x  G  B|c|  (5.2e) 


In  this  formulation  5.2a  is  the  objective  function,  5.2b  are  the  transshipment 
portion  of  the  flow  balance  constraints,  5.2c  are  the  source  and  sink  portion  of  the 
flow  balance  constraints,  5. 2d  are  the  knapsack  constraints,  and  5.2e  are  the  flow 
forcing  constraints. 

The  RCSP  is  a  variation  on  the  classic  Shortest  Path  Problem  (SPP)  in  which 
the  goal  is  to  find  a  path  of  minimum  total  weight  connecting  a  source  node,  s,  and 
a  sink  node,  t.  The  RCSP  adds  an  additional  knapsack  constraint  which  is  applied 
to  the  sum  of  the  weights  on  the  minimum  s,t-path.  Although  the  SPP  is  well 
solved  by  algorithms  like  Dijkstra’s  algorithm,  Floyd- Warshall,  and  the  Out-of-Kilter 
algorithm,  no  efficient  algorithm  for  solving  the  RCSP  has  yet  been  developed  [48]. 

Formulating  the  RCCP  in  multimodal  freight  networks  as  an  instance  of  the 
RCSP  was  first  accomplished  by  Hartlage  and  Weir  in  2012  [56].  The  original  for¬ 
mulation  is  shown  below  for  reference. 

We  first  define  the  following  sets  and  parameters  to  generalize  the  single  mode 
for  addressing  multimodal  freight  transportation.  Let  P  —  {1,  2, . . .  ,p}  be  the  set 
of  available  transportation  modes  and  k  E  P  be  be  the  index  for  these  modes.  Also 
N  =  {1,2,...,  77.}  be  the  set  of  nodes  representing  nodes,  as  in  the  single-mode 
formulation  with  i  E  N  being  the  index  for  the  nodes.  Expanding  the  original 
graph  as  described  in  the  previous  paragraph  results  in  a  node  set  with  np  nodes. 
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The  indices  for  the  node  set  of  the  expanded  graph  are  l,v  G  {1  Let 

Nk  =  {1  +  n(k  —  1 ),...,  n  +  n(k  —  1)}  be  the  node  set  for  mode  k  in  the  expanded 
graph.  Finally,  let  Ml  =  {i, . . .  ,i  +  n(k  —  1), . . . , i  +  n(p  —  1)}  be  the  node  set  in  the 
expanded  graph  corresponding  to  node  i  in  the  original  graph.  Define  the  additional 
sets,  and  matrices  used  in  formulating  the  multimodal  RCSP: 


Ck  =  {( l,v )  :  ( l,v )  is  an  intra-mode  arc  in  the  graph  for  mode  k} 

Dl  =  {( l,v )  :  ( l,v )  is  an  inter-mode  arc  in  the  graph  for  node  i} 
i]  =  (J  Nk\  the  node  set  for  the  expanded  graph 

k 

A  =  ^|J  Ck^j  1J  ^|J  D^j  :  the  arc  set  for  the  expanded  graph  (5-3) 

Gk  =  ( Nk,Ck ):  the  graph  for  mode  k 
H1  =  ( Ml ,  Dl):  the  graph  for  node  i 
T  =  (77,  A):  the  expanded  graph 

Now  we  define  additional  indices.  Let  u  G  {1, ...,  |A|}  be  the  index  for  the  arc 
set  of  the  expanded  graph,  mk  G  {1,  ...\Cl  |}  be  the  index  of  the  intra-mode  arc  set 
for  for  each  transportation  mode  k ,  and  zl  G  {1, ...,  |Dl|}  be  the  index  of  the  inter¬ 
mode  arc  set  for  each  node  i.  Also  let  s  and  t  be  the  indices  of  the  source  and  sink 
nodes  in  the  original,  non-expanded  graph  respectively.  With  these  sets  and  indices 
we  can  now  develop  the  multimodal  formulation  using  the  matrices  defined  below. 


80 


1 

if  /  G  Nk  is  tail  of  arc  mk  e  Ck 

rpk  _ 

< 

-1 

if  /  G  Nk  is  head  of  arc  mk  G  C 

0 

otherwise 

1 

if  /  G  Ml  is  tail  of  arc  zl  G  D' 

F*  — 

r  l,zx 

< 

-1 

if  /  G  Ml  is  head  of  arc  zl  G  Dl 

0 

V 

otherwise 

E  = 

[E1:.E2: . 

. . : Ep } 

F  = 

[EV:E2:. 

. .  \Fn] 

A  = 

[e:f] 

qs(u)  =  Y  Ai,uVu 

ieMs 
Ai,u  =  1 

q-t(u)  =  Y  A  i,uVu 

leM* 

Al,u=~l 


Q  = 


q-t 


(5.4) 


Additionally,  we  define  the  following  vectors  and  matrices.  Note  that  in  the 
definitions  below,  the  notation  d!  is  taken  to  represent  the  matrix  transpose  of  d. 


81 


h  :  the  1  x  p  matrix  of  mode  transportation  costs  in  dollars  per  MTM 

r  :  the  1  x  p  vector  of  MTM  available  in  mode  k 

B_sk  :  the  1  x  \Ck\  vector  of  arc  weights  for  mode  k 


B  tl  :  the  1  x  \D'l\  vector  of  arc  weights  for  node  i 

B_s  =  [B  s^-.B  s2\ . . .  \B_sp] 

B_t  =  [B_tr:B_t2\  .  .  .  \B_tn } 

x  :  the  1  x  |A|  vector  of  decision  variables 

\B  s'] 


B'  = 

B_t' 


J  =  [JV.  0] 


Using  the  preceding  developments  we  formulate  the  multimodal  RCSP  as: 


Minimize  Bx 
subject  to:  Ax  =  0 

Qx  =  1  (5.6) 

Jx  <  r 
x  e  B|A| 

This  multimodal  formulation  of  the  RCSP  is  the  focus  of  the  methodology  that 
follows.  Similar  research  efforts  have  used  metaheuristics  to  solve  a  similar  problem 
called  the  resource  constrained  project  scheduling  problem  (RCPSP).  The  RCPSP 
is  similar  to  the  RCSP  in  that  it  includes  constraints  on  total  available  resources 
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but  differs  in  that  its  objective  is  to  minimize  the  total  time  required  to  complete  a 
series  of  tasks  (called  the  “makespan”)  and  also  includes  precedence  constraints  on 
the  tasks. 

Calhoun  used  a  modified  version  of  Van  Hove’s  original  Multi-modal  resource 
constrained  project  scheduling  problem  with  generalized  precedence  constraints  (MMG- 
PRCPSP)  to  model  the  assignment  of  aircraft  to  targets.  He  formulated  the  MMG- 
PRCPSP  as  a  lexicographic  ordered  goal  program  to  minimize  the  number  of  targets 
not  covered.  He  solved  the  problem  using  Tabu  Search  [25,62],  In  the  context  of 
Calhoun  and  Van  Hove’s  research  “multi-modal”  refers  to  the  number  of  different 
ways  or  “modes”  in  which  a  particular  task  can  be  accomplished.  This  is  an  im¬ 
portant  distinction  since  the  term  “multi-modal”  is  used  in  this  article  to  mean  the 
movement  of  freight  using  two  or  more  modes  of  transportation. 

Cocllo  Coello  et  al.  develop  a  method  for  incorporating  multiple  objectives 
into  particle  swarm  optimization  (PSO)  by  generating  uniformly  distributed  points 
along  a  Pareto  front  [27]. 

Nasiri  develops  a  pseudo  PSO  metaheuristic  for  solving  the  RCPSP  [94],  The 
metaheuristic  uses  path  relinking  to  allow  particles  to  move  toward  optima  (local 
and  global). 

Merkle  et  al.  developed  an  Ant  Colony  Optimization  metaheruistic  for  the 
RCPSP  [84], 

5.2  Methodology 

The  methodology  developed  here  is  intended  to  support  tradeoff  analysis  in 
capacity  planning  decisions  as  they  relate  to  fleet  mix,  total  cost,  and  time.  In 
particular,  the  first-of-kind  methodology  developed  here  allows  for  variable  time  res¬ 
olution  modeling  (selection  of  time  interval  length  and  units),  and  provides  insight 
into  the  value,  in  dollars,  of  additional  time  or  transportation  assets.  Typical  ques¬ 
tions  that  may  be  addressed  by  the  methodology  developed  here  are: 
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•  What  is  the  minimum  number  of  days,  weeks,  or  months  in  which  to  ship  x 
tons  of  freight  with  the  specified  mix  of  assets? 

•  What  is  the  lexicographic  minimum  asset  mix  required  to  ship  x  tons  in  y  days, 
weeks,  or  months? 

•  What  is  the  minimum  cost  asset  mix  to  ship  in  a  given  number  of  days,  weeks, 
or  months? 

•  Is  the  current  asset  mix  minimal?  That  is,  are  there  excess  assets  in  the  current 
asset  mix  that  could  be  re-purposed  for  another  mission? 

Although  there  are  certainly  other  questions  that  can  be  addressed  using  this 
new  methodology,  these  typify  questions  a  decision  maker  might  be  interested  in 
regardless  of  whether  the  application  is  to  military  or  civilian  transportation. 

5.2.1  Finding  Near- Efficient  Solutions.  In  theory  the  Pareto  front  is  con¬ 
tinuous  and  therefore  contains  an  infinite  number  of  efficient  solutions.  In  some 
applications,  such  as  the  application  to  multimodal  transportation  discussed  in  this 
paper,  the  Pareto  front  is  discrete.  However,  the  number  of  discrete  Pareto  points 
on  the  frontier  can  be  quite  large.  Consider  a  relatively  small  number  of  transporta¬ 
tion  modes  and  number  of  assets  in  each  mode.  Given  7  transportation  modes  and 
20  of  each  transportation  asset  type,  the  number  of  potential  efficient  solutions  is 
20'  =  1,280,000,000.  Finding,  and  comparing  all  of  these  points  is  time  consum¬ 
ing  and  for  large  enough  problems  quickly  becomes  impractical.  Typically,  there  are 
three  items  to  consider  when  evaluating  multi  objective  optimization  problems  [118]. 

1.  Minimization  of  the  distance  between  the  near-Pareto  points  found  and  the 
global  Pareto  front. 

2.  Maximize  the  spread  of  the  near-efficient  solutions  found  to  provide  a  smooth 
and  uniform  near-Pareto  front. 

3.  Maximize  the  number  of  elements  of  the  global  Pareto  front  found. 
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We  use  an  ACS  metaheuristic  developed  by  Hartlage  and  Weir  to  aid  in  gener¬ 
ating  solutions  for  this  article  [55] .  The  results  included  an  optimization  of  the  ACS 
input  parameters  to  jointly  optimize  both  the  running  time  and  the  solution  quality. 
These  previous  results  suffice  for  addressing  the  first  and  third  item  in  the  list  above. 
In  this  article,  our  focus  is  on  finding  uniformly  distributed,  near-efficient  solutions. 

5.2.2  Distance  between  Solutions  and  Coverage  of  the  Pareto  Front.  By 
way  of  an  example,  consider  route  choices  generated  by  an  internet-based  route 
mapping  engine.  Initially,  a  start  and  an  end  point  are  specified.  The  engine  provides 
a  driving  route  between  the  specified  origin  and  destination.  However,  alternative 
routes  are  often  specified.  In  order  to  ensure  that  an  alternate  route  differs  from  the 
original  route  by  more  than  simply  cutting  through  a  parking  lot,  the  engine  must 
employ  a  measure  of  distance  between  successive  solutions. 

This  article  considers  a  two  dimensional  objective  space.  The  objectives  con¬ 
sidered  are  time  (in  days)  and  cost  (in  dollars).  For  a  fixed  quantity  of  freight,  there 
is  an  inherent  tradeoff  between  delivery  time  and  cost.  Increasing  time  reduces  the 
number  of  transportation  assets  and  the  associated  cost  of  capacity.  Decreasing  the 
time  requires  the  use  of  more  transportation  assets  and  is  likely  to  involve  faster 
and  more  expensive  transportation  modes,  measured  in  dollars  per  million-ton-mile, 
to  meet  shorter  timelines.  Shipping  freight  by  air  is  an  example  of  a  relatively  fast 
but  expensive  mode  of  transportation.  Increasing  time  will  decrease  the  number  of 
transportation  assets  required  to  provide  required  capacity  and  will  tend  to  require 
less  expensive  modes  of  transportation. 

We  can  measure  how  uniformly  the  solutions  are  distributed  in  the  objective 
space  by  measuring  the  range  variance  of  the  near-efficient  solutions  in  the  objective 
space.  The  metric  for  measuring  range  distance  proposed  by  Schott  is  [103]: 
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n 


(5.7) 


S 


\ 


n  —  1 


£  (a  -  *)2 


Where  i  and  j  are  indices  of  the  total  number  of  efficient  solutions  generated, 
n,  and  di  is  the  distance  from  a  solution,  i,  to  the  nearest  solution  j  for  j  =  1,  ...,n. 
Mathematically,  d*  is  expressed  as: 


dj  =  mini 


+  f-2  ( )  —  fi  ( ) 


Vi,J  =  1,  2, ...,  n 


(5.8) 


Also,  let  d  be  the  arithmetic  average  of  all  d*  and  fk  be  vectors  for  k  =  1,  2, ...,  n. 
This  metric  allows  us  to  evaluate  the  uniformity  of  near-efficient  solutions.  Coverage 
of  the  Pareto  front  is  another  issue  which  requires  discussion.  Coverage  can  be 
determined  by  selecting  two  values  in  the  objective  space.  For  a  fixed  quantity  of 
freight,  selecting  minimum  and  maximum  values  for  time  determines  the  size  and 
location  of  the  Pareto  front  segment  of  interest.  Setting  a  minimum  value  for  time 
(Tmiri),  we  obtain  the  highest  cost  near-efficient  solution  (Cmax).  Setting  a  maximum 
allowable  time  ( Tmax ),  we  obtain  the  lowest  cost  near-efficient  solution  ( Cmin )• 

5.2.3  Generalized  Range  Variance  for  Lexicographic  Objective  Ordering. 
The  definition  of  range  variance  given  by  Schott  is  only  sufficient  to  measure  the 
range  variance  in  a  two  dimensional  objective  space.  In  the  case  of  lexicographic 
ordering  of  m  components  we  require  range  variance  measured  in  all  m  dimensions. 
Generalizing  range  distance  requires  us  to  redefine  the  quantity  d*  as: 


di  =  minj 


Vi,  j 


1,2, ...,  n 


(5.9) 


Defining  d*  in  this  way  we  can  now  provide  range  variance  not  only  in  in¬ 
stances  where  the  objective  space  has  the  dimensions  of  cost  and  time  but  also  in 


situation  where  we  are  interested  in  an  objective  space  having  dimensions  of  time 
and  lexicographic  order  of  all  transportation  asset  types. 

5.2.4  The  Value  of  Dominated  Solutions.  In  the  context  of  capacity  plan¬ 
ning  for  multimodal  transportation,  global  Pareto  solutions  have  no  excess  capacity. 
Stated  in  another  way,  the  loss  of  a  single  asset  from  among  the  fleet  proposed  in  the 
Pareto  solution  causes  the  problem  to  be  come  infeasible  for  the  current  value  of  T . 
In  some  applications  and  situations,  a  decision  maker  may  wish  to  include  reserve 
capacity  in  the  planning  phase  in  order  to  account  for  reliability  issues  with  the  fleet 
or  late  departures  of  assets  due  to  weather  delays  for  example.  Making  a  decision 
to  operate  with  excess  capacity  is  equivalent  to  selecting  a  dominated  point  (non¬ 
efficient  solutions)  in  the  objective  space.  Therefore,  presentation  of  non-efficient 
solutions  is  beneficial  to  decision  makers  desiring  this  additional  capacity. 

5.2.5  Methodology  Development.  The  input  required  for  this  methodology 
include  parameters  and  data  which  describe  the  transportation  mode  characteristics 
and  the  configuration  of  the  underlying  network.  Also,  the  maximum  number  of  each 
transportation  asset  type  (mode)  is  a  required  input.  An  overview  of  the  process  is 
provided  here  for  reference. 

•  Provide  required  parameter  and  data  inputs. 

•  Determine  the  minimum  number  of  days  required  to  ship  the  required  tons  of 
freight  with  the  provided  asset  mix. 

•  Using  data  provided  for  number  and  length  of  time  intervals,  generate  near¬ 
efficient  solutions  for  each  of  the  interval  endpoints. 

•  Present  the  full  set  of  near-efficient  solutions  to  the  decision  maker. 
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5.2.5. 1  Parameters  and  Data.  There  are  a  number  of  inputs  that  are 
required  to  initialize  the  problem.  The  following  parameters  determine  characteris¬ 
tics  of  the  network. 

•  hk'-  The  transportation  cost,  in  dollars,  per  million-ton-mile  for  mode  k. 

•  LexObjOrder:  The  lexicographic  objective  ordering  for  all  modes  of  trans¬ 
portation  [1  :  m\. 

•  UTErate :  The  fleet  utilization  rate,  in  hours,  for  each  transportation  mode. 

•  Productivity:  The  fleet  productivity,  or  percentage  of  the  time  that  assets  do 
not  move  while  empty  (for  repositioning)  for  each  transportation  mode. 

•  BlockSpeed:  The  average  speed,  door-to-door,  for  each  transportation  mode. 

•  Payload:  The  single  trip,  payload  of  each  transportation  mode  in  tons. 

•  XFerCost:  The  cost,  per  million  tons,  to  transfer  freight  between  two  modes 
of  transportation  (could  be  unique  for  each  mode  pairing). 

•  s:  The  source  node  or  point  of  embarkation  (POE)  for  the  system. 

•  t:  The  sink  node  or  point  of  debarkation  (POD)  for  the  system. 

•  rc:  The  required  system  capacity  in  millions  of  tons. 

•  IntSize:  The  number  of  days  in  a  time  interval. 

•  IntNum:  The  number  of  time  intervals  to  analyze. 

•  Max  Assets:  A  vector,  [cq,  02, ...,  am]  of  the  number  of  transportation  assets 
available  in  each  mode  of  transportation. 

•  T:  A  scalar  denoting  the  number  of  days  allowed  to  ship  the  amount  of  freight 
specified  by  rc. 

The  next  set  of  parameters  determine  the  search  characteristics  for  the  ACS- 


RCSP  metaheuristic: 


•  N umN eighbor s:  The  number  of  nearest  neighbors  on  the  reduced  nearest 
neighbor  (NN)  list.  A  reduced  NN  list  reduces  computation  time  by  as  much 
as  ten  percent  [37]. 

•  a:  pheromone  weight  factor  on  the  interval  [0, 1]  reduces  influence  of  the 
pheromones. 

•  (3:  heuristic  weight  factor  in  the  interval  [0,1],  a  lower  (3  reduces  the  influence 
of  the  heuristic  information. 

•  To  weight  parameter  for  pheromone  update  operation 

•  g0:  An  intensification/diversification  parameter  in  the  interval  [0,1],  a  higher 
go  means  more  exploitation  of  learned  knowledge  while  a  lower  g0  means  more 
exploration  of  arcs. 

5. 2. 5. 2  Parameters  and  Transportation  Modes.  We  can  conveniently 
represent  the  time  and  transportation  asset  input  variables  for  each  problem  instance 
as  a  single  vector.  By  convention,  the  first  variable  in  the  vector  is  the  parameter  T, 
in  days.  The  remaining  values  in  the  vector  indicate  how  many  of  each  transportation 
mode  asset  type  is  available.  These  variables  are  captured  in  Max  Assets  defined 
in  the  previous  section.  For  the  purposes  of  this  methodology,  a  transportation 
mode  is  a  grouping  of  homogeneous  transportation  assets  with  access  to  a  common 
network  structure  (nodes  with  adjacencies  defined).  Defining  modes  in  this  way 
allows  greater  flexibility  in  accurately  describing  the  transportation  system.  For 
instance,  although  it  might  be  tempting  to  classify  all  air  transportation  assets  as 
belonging  to  a  single  mode  of  transportation,  such  a  classification  does  not  accurately 
capture  the  differences  between  the  asset  types 

5. 2. 5. 3  Determining  Minimum  Days  Required.  Determining  the  min¬ 
imum  number  of  days  required  to  provide  necessary  system  capacity  involves  produc¬ 
ing  an  infeasible  solution  and  using  bisection  to  measure  the  infeasibility  in  terms  of 


days.  Producing  an  infeasible  solution  is  relatively  easy.  However,  we  are  interested 
in  generating  an  infeasible  solution  using  the  maximum  number  of  transportation 
assets  available  in  all  modes  of  transportation.  To  do  this,  we  set  T  =  1  and  solve  the 
problem.  From  this  infeasible  solution,  we  implement  a  reverse  and  forward  bisection 
routine  to  determine  the  minimum  number  of  days,  T,  required  to  obtain  a  feasible 
solution  with  Max  Assets.  First,  we  must  define  a  function  which  determines  the 
feasibility  of  a  RCSP  provided  inputs  of  T  and  Max  Assets  as  follows: 


/( Time ,  Asset  Mix) 


1  if  RCSP  admits  a  feasible  solution 

—  1  otherwise 


(5.10) 


The  function  /  used  in  this  methodology  is  ACS-RCSP  or  a  binary  integer 
programming  solver.  ACS-RCSP  is  the  Ant  Colony  System  metaheuristic  for  solving 
the  RCSP  originally  developed  by  Hartlage  and  Weir  [55].  The  reverse  and  forward 
bisection  routine  pseudo-code  is  described  in  Figure  1: 

This  algorithm,  FBintBisect,  is  similar  to  the  bisection  method  used  to  find 
roots  of  a  function  [24],  however,  it  has  been  modified  to  include  two  additional 
features.  First,  it  includes  a  preprocessing  step  that  “reverse”  bisects,  increasing 
the  value  of  T  until  a  feasible  solution  is  returned.  Second,  it  uses  the  function 
ceiling(x )  which  rounds  fractional  values,  x ,  to  the  nearest  integer  value  in  the 
direction  of  positive  infinity. 

5. 2. 5. 4  Generating  Near- Efficient  Solutions.  The  value  returned  by 
the  bisection  algorithm,  minT.  is  the  minimum  number  of  days  required  in  order 
for  the  problem  to  be  feasible  with  the  assets  available  in  MaxAssets.  Now,  be¬ 
ginning  from  minT  we  generate  the  number  of  near-efficient  solutions  specified  by 
IntSize  and  IntNum.  The  subroutine  to  generate  these  near-efficient  solutions 
requires  first  developing  three  other  utility  subroutines:  “StVallnc” ,  uAssetMixn , 
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Algorithm  1  Forward-Backward  Integer  Bisection  Subroutine 
L  <r~  0 
H 

x  f(H,  Assets ) 

while  x  =  —  1  do 
H  ^2*H 
x  f(H,  Assets ) 

end  while 

M  ceiling {{L  +  H)/ 2) 

while  i/  —  L  >  2  do 
if  a;  =  1  then 

H  4r-  M 

M  ceiling((L  +  H)/ 2) 

else  if  a:  =  —  1  then 
L  M 

M  •(—  ceiling ((L  +  H)/2) 

end  if 

x  •(—  f(H ,  Max  Assets) 

T  <r-  M 

end  while 

return  T 


and  “Feas  Asset  Min" .  StVallnc  takes  an  infeasible  solution  and  determines  by  how 
much  each  of  the  transportation  asset  types  and  time  would  have  to  be  increased 
independent  of  one  another  in  order  to  achieve  feasibility.  In  essence,  StVallnc  pro¬ 
vides  a  measure  of  by  how  much  the  problem  is  infeasible  with  respect  to  each  of  the 
controllable  inputs.  AssetMix  increases  system  capacity  by  adding  assets  (up  to  that 
specified  in  MaxAssets)  in  non-decreasing  order  with  respect  to  the  lexicographic 
ordering  (or  by  cost  per  million  ton  mile  if  no  objective  ordering  is  provided)  until 
a  feasible  solution  is  returned.  Once  this  feasible  solution  is  found,  Feas Asset M in 
determines  if  the  asset  mix  returned  by  AssetMix  contains  any  excess  assets.  Any 
excess  assets  are  removed  from  the  asset  mix  so  that  the  final  asset  mix  returned  by 
Feas  As  set  Ad in  is  minimal,  in  that  it  is  both  necessary  and  sufficient  for  feasibility 
at  the  current  value  of  T.  StVallnc  is  described  in  Algorithm  2,  FeasAssetMin  in 
Algorithm  4,  and  AssetMix  in  Algorithm  3. 
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Algorithm  2  Subroutine  for  Determining  Initial  Increase  to  Attain  Feasibility 

StVallnc  zeros{  1,  length(MinAssets)) 

for  i  —  1  :  m  do 

StValInc{i)  FBintBisect(Time,  Min Assets) 

end  for 

return  StVallnc _ 

For  convenience  in  expressing  As  set  M  ix  in  pseudocode  define  the  following 
variables: 

•  R(x ):  The  number  of  additional  assets  in  mode  x  required  to  make  the  problem 

feasible. 

•  A(x):  The  number  of  additional  assets  in  mode  x  available. 

•  U(x):  The  number  of  assets  in  mode  x  currently  used  (Note:  total  assets  in 

mode  x  =  A(x)  +  U(x)). 

Now,  let  AdinAssets  be  an  m-vector  of  ones.  Notice  that  the  choice  of  IntSize 
and  IntNum  determine  the  distance  between  efficient  solutions  and  the  number  of 
solutions  generated,  respectively.  When  combined,  the  selection  of  these  parame¬ 
ters  determines  both  fidelity  and  coverage  of  the  decision  space.  The  range  variance 
provides  a  measure  of  the  uniformity  of  distribution  of  the  points  generated  in  the  ob¬ 
jective  space.  The  routine  for  generating  efficient  solutions  is  described  in  Algorithm 
5. 

5.3  Example  Problems 

We  demonstrate  the  method  of  finding  near-efficient  solutions  in  this  section. 
Results  are  presented  as  a  continuum  from  the  most  costly,  shortest  time  required 
solution  to  the  least  costly  most  time  required  solution.  For  each  solution,  the  path, 
cost  in  dollars,  and  the  asset  mix  are  provided. 

5.3.1  Randomly  Generated  Examples.  The  test  network  consists  of  50 
nodes,  three  modes  of  transportation,  and  the  node  positions  are  randomly  generated 
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Algorithm  3  Subroutine  for  Determining  the  Asset  Mix 
if  lexOrd  =  1  then 

while  “no  feasible  solution  is  returned”  do 

x  =  “index  of  least  important  transportation  mode  such  that”  A(x)  >  0 
and  R(x)  >  0 

if  R(x)  <=  A(x)  and  A(x)  >  0  then 
U{x)  <-  U(x)  +  R(x) 

A(x)  <—  A(x)  —  R(x) 
else  if  R(x)  >  A(x)  and  A(x)  >  0  then 
U(x)  <-  U(x)  +  A(x) 

A{x)  0 

end  if 

for  i  —  1  :  m  do 

R(i)  =  StValInc(Time,U(i )) 

end  for 
end  while 

else  if  lexOrd  =  0  then 

while  “no  feasible  solution  is  returned”  do 

x  =  “index  of  least  expensive  transportation  mode  such  that”  A{x)  >  0 
and  R(x)  >  0 

if  R(x)  <=  A(x)  and  A(x)  >  0  then 
U(x)  <-  U(x)  +  R(x) 

A(x)  A(x)  —  R(x) 
else  if  R(x)  >  A(x)  and  A(x)  >  0  then 
U(x)  <-  U(x)  +  A{x) 

A{x)  0 

end  if 

for  i  —  1  :  m  do 

R(i)  =  StValInc(Time,U(i )) 

end  for 
end  while 
end  if 

return  U  (x) 
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Algorithm  4  Subroutine  for  Determining  Minimal  Assets 
path  =  “sequence  of  nodes  on  the  feasible  s,t-path” 

L  =  “number  of  nodes  in  path 

I A  =  “the  vector  of  initial  assets  used  in  each  mode  1  :  m” 

EA  =  “the  vector  of  excess  assets  used  in  each  mode  1  :  m” 

RA  =  “the  vector  of  reduced  assets  used  in  each  mode  1  :  m" 

EM  =  “the  vector  of  excess  MTM  used  in  each  mode  1  :  m” 

AM  =  “the  vector  of  available  MTM  used  in  each  mode  1  :  m” 

IM  =  “the  vector  of  initial  MTM  used  in  each  mode  1  :  m” 

RM  =  “the  vector  of  reduced  MTM  used  in  each  mode  1  :  m” 

U  =  “the  vector  of  UTE  rates  for  each  mode  1  :  m” 

PA  =  “the  vector  of  single-trip  payloads  (in  tons)  for  each  mode  1  :  m” 
PP  =  “the  vector  of  productivity  rates  for  each  mode  1  :  m” 

BL  =  “the  vector  of  block  speeds  for  each  mode  1  :  m” 

AP  =  “the  vector  of  max  assets  available  for  each  mode  1  :  m” 
for  j'=l:  m  do 

for  i  =  1  :  L  —  1  do 

if  mode  of  path(i)  =  mode  of  path(i  +  1)  =  j  then 
IM(j )  IM(j)  +  arclength  (path(i) , path(i  +  1)) 

end  if 
end  for 
end  for 

for  %  =  1  :  m  do 

MTMasset{i)  (P (i)  *  PA(i)  *  BL(i)  *  PP(i))/l,  000,  000 
AM(i)  MTMasset(i )  *  AP(i)  *  Time 

end  for 

EM  <r-  AM  -  IM 

for  i  =  1  :  m  do 

EA(i)  =  floor(A™}2ime) 

end  for 

for  i  —  1  :  m  do 

RA(i)  =  IA(i)  -  EA(i) 

end  for 

return  RA 
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Algorithm  5  Subroutine  for  Generating  Near-Efficient  Solutions 
if  “lexicographic  ordering  is  provided”  then 
lexOrd  =  1 

else  if  “lexicographic  ordering  is  not  provided”  then 
lexOrd  =  0 

end  if 

for  i  —  1  :  IntNum  +  1  do 
R  =  minT  +  i  *  IntSize 
x  =  f(R,  Min  Assets) 

if  x  =  —  1  then 

minDays  =  StVallnc 

end  if 

runDays  =  minDays  +  intSize  *  (i.  —  1 ) 

AssetMixVect  =  AssetMix(lexOrd,  runDays) 
minAssetMixVect  =  FeasAssetMin(AssetMixVect) 

end  for 


on  a  one  thousand  square  mile  grid.  All  parameter  values  are  summarized  in  the 
Table  5. 3. 1.1: 


5.3. 1.1  An  Example  Problem.  The  results  of  the  analysis  are  provided 
in  tables  5.3  and  5.4  below.  Notice  that  the  minimum  number  of  days  required  to 
provide  the  required  capacity  with  the  given  assets  is  127.  This  value  is  calculated 
using  all  available  assets.  Another  feature  of  the  methodology  illustrated  by  the 
results  is  that  the  asset  mix  for  the  minimum  number  of  days  is  less  than  all  available 
assets.  In  fact,  using  bisection,  this  asset  mix  is  minimal.  In  other  words,  a  reduction 
of  a  single  unit  of  any  asset  in  the  given  asset  mix  results  in  an  infeasible  problem 
for  a  given  value  of  T. 

Two  other  concepts  bear  explanation.  First,  notice  that  costs  of  solutions  3, 
4,  and  5  in  the  asset/time  example  are  identical.  The  cost  expressed  is  for  the 
number  of  MTM  used  rather  than  the  number  available  for  the  given  asset  mix. 
The  increase  in  time  between  solution  one  and  solution  two  results  in  an  increase  in 
MTM  provided  on  a  per  asset  basis.  Therefore,  although  the  number  assets  in  the 
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interval  length  (in  days) 

30 

number  of  intervals 

5 

grid  size  (in  square  miles) 

800 

POE  (source  node) 

1 

POD  (sink  node) 

50 

number  of  millions  of  tons  to  move 

5 

arc  density 

90% 

max  assets  available 

[10,15,10] 

cost  per  MTM  (in  $) 

[100,110,120] 

lexicographic  ordering  (3  =  most  important) 

[3,2,1] 

UTE  rate  (hours/day) 

[15,15,15] 

productivity  (%) 

[70,70,70] 

block  speed  (mph) 

[515,366,579] 

payload  (tons) 

[85,22,135] 

transfer  costs  (in  $  per  MT) 

0.000001 

Table  5.2:  Small  Experiment  Parameters 


asset  mix  has  decreased,  the  same  number  of  MTM  are  being  utilized  to  provide  the 
required  capacity. 

Secondly,  notice  that  the  cost  increases  in  subsequent  solutions  in  the  assets/- 
tirne  example  while  they  are  reduced  in  the  cost /time  example.  This  is  explained 
by  the  lexicographic  order  when  using  the  lexicographic  order  rather  than  cost  as 
an  objective  component,  the  most  expensive  mode  of  transportation  is  given  an  ob¬ 
jective  component  ordering  as  the  least  important  or  critical  transportation  mode. 
Therefore,  the  methodology  seeks  to  minimize  assets  in  the  other  two  modes  of  trans¬ 
portation  first  which  results  in  a  greater  amount  of  capacity  being  provided  by  the 
most  expensive  (but  least  important)  transportation  mode. 

Using  the  generalized  range  variance  developed  in  section  5.2,  the  range  vari¬ 
ance  for  these  solutions  is  calculated  as  0.0601  indicating  that  the  near-efficient  points 
are  uniformly  spaced  in  the  objective  space. 

Now  if  we  solve  the  same  problem  based  upon  minimizing  cost  rather  than 
based  upon  a  lexicographic  ordering  of  the  transportation  modes  we  obtain  the 
solutions  in  table  5.4. 
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solution 

days 

assets 

$  millions 

path  (original) 

path  (with  modes) 

1 

127 

[9,0,10] 

177.99 

[1,10,32,50] 

mode  1  mode  3  model 

[1,10]  [10,  32]  [32,  50] 

2 

157 

[7,0,8] 

178.09 

[1,32,50] 

mode  3  mode  1 

[1,32]  [32,  50] 

3 

187 

[5,0,10] 

187.95 

[1,10,50] 

mode  1  mode  3 

[1,10]  [10,  50] 

4 

217 

[5,0,9] 

187.95 

[1,10,50] 

model  mode  3 

[1,10]  [10,  50] 

5 

247 

[4,0,8] 

187.95 

[1,10,50] 

model  mode  3 

[1,10]  [10,  50] 

6 

277 

[0,0,7] 

188.03 

[1,50] 

mode  3 

[1,50] 

Table  5.3:  Multiobjective  Lexicographic  Order  &  Time  Solution  Results 


solution 

days 

assets 

$  millions 

path  (original) 

path  (with  modes) 

1 

127 

[9,0,10] 

177.99 

[1,10,32,50] 

mode  1  mode  3 

[1,10]  [10,32,50] 

2 

157 

[10,0,7] 

174.53 

[1,44,50] 

model  mode  3 

[1,44]  [44,  50] 

3 

187 

[9,0,4] 

170.76 

[1,44,50] 

model  mode  3 

[1,44]  [44,  50] 

4 

217 

[9,0,2] 

166.69 

[1,44,50] 

model  mode  3 

[1,32]  [32,  50] 

5 

247 

[9,0,2] 

164.60 

[1,44,32,50] 

model  mode  3  model 

[1,44]  [44,  32]  [32,  50] 

6 

277 

[9,0,2] 

163.22 

[1,44,32,50] 

model  mode  3  model 

[1,44]  [44,  32]  [32,  50] 

Table  5.4:  Multiobjective  Cost  &  Time  Solution  Results 
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The  range  variance  of  this  solution  is  calculated  as  1613.09.  This  range  vari¬ 
ance  indicates  less  uniformity  in  the  spacing  of  the  efficient  solutions.  However,  the 
magnitude  of  the  cost  objective  component  is  quite  large  numbering  in  the  hundreds 
of  thousands  whereas  the  magnitude  of  the  largest  objective  component  of  time  in 
the  cost /time  example  was  only  in  the  hundreds. 

Each  node  in  the  network  is  given  a  unique  identifier  which  indicates  the  cor¬ 
responding  node  and  mode  number.  Each  node  v  —  %  +  n{k  —  1)  where  i  is  the  node 
identifier  in  the  original  network  and  k  is  the  number  of  modes  of  transportation 
in  the  problem  [56].  For  example,  in  a  multimodal  network  with  n  =  10  nodes  and 
k  =  3  modes  of  transportation,  node  five  in  mode  three  would  receive  the  unique 
identifier  v  —  5  +  10(3  —  1)  =  25.  With  this  node  labeling  method  we  describe  one 
of  the  multimodal  paths  in  the  solutions  above.  The  remaining  paths  can  easily  be 
obtained  by  the  reader.  The  path  in  solution  1  in  Table  5.3  consists  of  three  shipping 
arcs  and  two  mode  transfer  arcs.  Arc  (1, 10)  is  in  mode  1.  At  node  10  capacity  is 
transferred  to  mode  3.  Capacity  is  added  to  arc  (10,  32)  in  mode  3.  At  node  32 
capacity  is  transfered  to  mode  1  and  capacity  is  added  to  arc  (32,  50)  completing  the 
path. 

5.4  Conclusions  and  Future  Work 

5. 4-1  Conclusions.  The  methodology  developed  in  this  research  provides 
insight  into  the  interaction  of  time,  transportation  asset  quantity  requirements,  and 
total  cost.  Since  the  RCSP  is  formulated  based  on  a  desire  to  generate  a  rough-cut 
capacity  estimate,  the  data  inputs  required  to  complete  this  analysis  are  limited  to 
data  that  are  almost  always  available  and  will  not  need  to  be  collected  such  as  asset 
capacity,  speed,  utilization  rate,  and  productivity  among  others. 

Very  little  input  is  required  on  the  part  of  the  decision  maker  in  order  to  com¬ 
plete  the  analysis.  Once  the  distribution  network  is  captured,  and  the  engineering 
data  for  all  transportation  assets  is  obtained,  analysis  only  requires  three  pieces  of 


information  from  a  decision  maker:  maximum  number  of  transportation  assets  avail¬ 
able  in  each  mode,  length  of  an  interval  (in  days)  between  solutions,  and  number  of 
intervals  for  which  to  generate  solutions.  The  latter  two  items  determine  both  the 
fidelity  and  solution  space  coverage. 

This  methodology  provides  an  explicit  analysis  of  the  tradeoffs  involved  be¬ 
tween  cost,  time,  and  transportation  assets  and  provides  a  decision  maker  with  the 
capability  to  determine  a  point  on  this  continuum  at  which  an  acceptable  balance 
of  all  three  goals  is  achieved. 
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VI.  Summary  &  Conclusions 

This  dissertation  has  made  several  original  contributions  to  the  held  of  Operations 
Research.  First,  a  new  methodology  for  modeling  the  multimodal  RCCP  using  the 
RCSP  was  developed  in  Chapter  3  along  with  results  showing  the  validity  of  the 
method  using  a  binary  integer  programming  formulation  and  a  commercial  solver. 

Chapter  4  developed  a  new  ACS  metaheuristic  for  solving  the  RCSP  which  is 
a  NP-Complete  optimization  problem.  Results  of  this  section  show  that  the  meta¬ 
heuristic  is  able  to  quickly  find  high  quality  feasible  solutions  to  large  instances  of 
the  RCSP. 

Chapter  5  builds  upon  the  modeling  approach  of  Chapter  3  and  the  ACS 
metaheuristic  of  chapter  4  by  developing  a  new  metaheuristic  for  generating  near 
Pareto  solutions  to  the  multimodal  RCCP.  The  methodology  considers  the  situation 
where  cost  and  time  are  the  objective  components  of  interest  and  also  considers 
a  lexicographic  objective  ordering  by  which  modes  of  transportation  can  be  rank 
ordered  according  to  importance  or  criticality.  An  extension  of  the  range  variance 
equation  is  developed  to  analyze  the  uniformity  of  near  efficient  solution  spacing  on 
the  near  Pareto  front  generated. 

Future  research  should  focus  on  generating  larger  instances  of  the  RCSP  for 
formulation  as  a  BIP  and  solution  by  commercial  solver  software.  The  focus  of  this 
dissertation  is  on  rough-cut  capacity  planning  which  seeks  to  answer  the  question: 
“Is  there  sufficient  capacity  in  the  system  to  meet  transportation  demands.”  An¬ 
other  area  of  future  research  is  aggregate  transportation  planning.  In  production 
and  inventory  management  this  planning  concept  is  used  to  answer  the  question: 
“How  much  should  I  produce  in  each  planning  period  in  order  to  meet  demands 
across  all  planning  periods.  Aggregate  transportation  planning  would  ask  an  analo¬ 
gous  question:  “What  quantities  of  freight  should  be  transported  in  each  period  to 
meet  demands  across  all  planning  periods.”  Aggregate  and  rough-cut  capacity  plan- 
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ning  are  closely  related  since  the  system  capacity  presents  a  constraint  that  must  be 
considered  in  aggregate  planning.  Aggregate  planning  in  transportation  should  con¬ 
sider  crew  scheduling,  labor,  warehousing,  inventory,  production,  shipping,  seasonal 
demand,  etc. 

The  author  believes  this  area  of  research  to  be  particularly  rich  it  topic  areas, 
and  timely  considering  the  current  need  to  reduce  expenditures  while  maintaining 
levels  of  transportation  capacity.  One  way  to  achieve  this  is  through  the  application 
and  extension  of  current  Operations  Research  methods  in  order  to  improve  efficiency 
and  effectiveness  of  transportation  planning  systems. 
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VII.  Appendix  1:  FCD  Experiment  Results  For  ACS  Testing 


102 


alpha 

beta 

maxjterations 

num_ants 

nummeighbors 

q_0 

rho 

xi 

percent  over  optimal 

run  time 

-1 

1 

-1 

1 

-1 

-1 

-1 

1 

0.820803265 

13.89314382 

-1 

1 

1 

-1 

1 

-1 

1 

1 

0.004476521 

11.47112901 

1 

-1 

-1 

1 

1 

-1 

-1 

-1 

8.356031737 

10.39651338 

-1 

-1 

1 

-1 

1 

-1 

-1 

1 

0.004476511 

13.6202956 

-1 

-1 

-1 

-1 

-1 

-1 

-1 

1 

0.220307227 

1.160888884 

1 

1 

-1 

-1 

-1 

1 

-1 

-1 

0.004476511 

1.182832176 

1 

1 

-1 

1 

-1 

-1 

-1 

-1 

3.323575033 

8.186534322 

-1 

-1 

1 

-1 

-1 

1 

-1 

1 

0 

12.6164134 

-1 

-1 

1 

1 

-1 

1 

1 

1 

0 

284.0650584 

-1 

1 

-1 

-1 

1 

-1 

1 

-1 

0.004476511 

1.924434949 

0 

0 

0 

-1 

0 

0 

0 

0 

0.004476506 

3.033284957 

0 

0 

0 

0 

0 

0 

0 

0 

0.004476511 

7.532166937 

-1 

1 

1 

1 

1 

1 

-1 

-1 

0 

351.2718527 

-1 

1 

1 

-1 

1 

1 

1 

-1 

0 

13.5170625 

1 

-1 

-1 

1 

-1 

1 

1 

1 

0 

12.07036425 

0 

0 

1 

0 

0 

0 

0 

0 

4.07792742 

59.80128759 

-1 

-1 

1 

1 

1 

-1 

1 

1 

12.75990173 

105.9970411 

1 

1 

-1 

1 

1 

-1 

-1 

1 

0.004476516 

11.03437797 

1 

-1 

-1 

1 

-1 

-1 

1 

-1 

4.077927415 

7.05146101 

-1 

-1 

-1 

1 

-1 

-1 

-1 

-1 

4.07792743 

8.662932021 

1 

-1 

1 

-1 

1 

-1 

1 

1 

11.92272051 

12.79150227 

1 

-1 

-1 

-1 

1 

1 

1 

1 

0.004476511 

1.158262433 

1 

-1 

1 

-1 

1 

1 

-1 

1 

0 

12.25825894 

1 

1 

1 

1 

-1 

-1 

-1 

1 

0.82080326 

296.1937443 

-1 

-1 

-1 

1 

1 

-1 

-1 

1 

0.831086705 

11.71672769 

1 

-1 

-1 

-1 

-1 

-1 

1 

1 

4.07792742 

1.040556272 

1 

1 

1 

-1 

-1 

1 

-1 

1 

0 

12.9834465 

1 

-1 

-1 

-1 

1 

1 

-1 

-1 

0.004476511 

1.172273152 

-1 

1 

-1 

1 

-1 

1 

-1 

-1 

0 

13.38079281 

-1 

-1 

1 

1 

-1 

-1 

-1 

1 

3.323575038 

364.9890439 

1 

0 

0 

0 

0 

0 

0 

0 

3.338162012 

9.002755658 

1 

1 

-1 

1 

1 

1 

-1 

-1 

0 

13.71408643 

1 

1 

-1 

-1 

1 

-1 

-1 

-1 

18.49638851 

1.772449152 

-1 

1 

-1 

1 

-1 

-1 

1 

-1 

0.004476506 

12.65459619 

1 

1 

1 

1 

-1 

-1 

1 

-1 

0.004476511 

287.9672545 

1 

1 

-1 

1 

1 

-1 

1 

-1 

29.95927477 

10.06518692 

1 

1 

1 

-1 

1 

1 

-1 

-1 

0 

13.88551646 

1 

-1 

1 

1 

1 

1 

-1 

-1 

0 

349.5637215 

-1 

-1 

1 

1 

-1 

-1 

1 

-1 

3.323575033 

390.7650016 

-1 

1 

1 

1 

-1 

-1 

1 

1 

3.323575038 

389.5671404 

0 

0 

0 

0 

1 

0 

0 

0 

4.07792742 

16.55031827 

-1 

1 

1 

1 

-1 

1 

-1 

1 

0 

279.400979 

-1 

1 

1 

1 

1 

-1 

-1 

1 

20.14362899 

113.3697083 

-1 

1 

-1 

1 

1 

-1 

1 

1 

9.696973944 

10.61463901 

-1 

1 

-1 

-1 

-1 

1 

1 

-1 

0.004476511 

1.092379199 

-1 

1 

-1 

-1 

1 

-1 

-1 

1 

29.95927476 

2.044553847 

-1 

1 

1 

-1 

-1 

-1 

1 

-1 

26.19590945 

11.50338954 

1 

1 

1 

-1 

-1 

-1 

-1 

-1 

0.82080326 

13.74440145 

1 

-1 

-1 

1 

1 

1 

1 

-1 

0 

13.28733065 

1 

1 

1 

1 

1 

-1 

1 

1 

29.95927476 

110.4218117 

1 

1 

1 

1 

-1 

1 

1 

1 

0 

286.3668183 

Table  7.1:  FCD  Experimental  Results  Part  1 
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alpha 

beta 

maxjterations 

nunrcants 

nummeighbors 

q_0 

rho 

xi 

percent  over  optimal 

run  time 

1 

1 

-1 

1 

-1 

-1 

1 

1 

0.004476511 

8.776136073 

-1 

1 

-1 

-1 

1 

1 

-1 

-1 

0.004476511 

1.188426881 

1 

-1 

1 

-1 

1 

1 

1 

-1 

0 

14.14519527 

1 

-1 

-1 

-1 

-1 

1 

-1 

1 

0.004476511 

1.168844644 

0 

0 

0 

0 

0 

0 

1 

0 

4.07792743 

8.870251945 

1 

-1 

1 

1 

-1 

1 

1 

-1 

0 

347.4683644 

-1 

-1 

1 

-1 

-1 

-1 

-1 

-1 

4.725338158 

7.689125853 

1 

-1 

1 

-1 

-1 

1 

-1 

-1 

0 

13.23266337 

-1 

-1 

-1 

-1 

1 

1 

-1 

1 

0.004476511 

1.118806301 

-1 

1 

1 

1 

1 

-1 

1 

-1 

3.323575038 

101.132139 

1 

1 

-1 

-1 

-1 

-1 

1 

-1 

4.077927425 

1.021030606 

1 

-1 

1 

1 

1 

-1 

-1 

1 

8.356031737 

90.26189041 

0 

0 

0 

0 

0 

-1 

0 

0 

28.9359434 

9.238635357 

1 

-1 

-1 

1 

1 

-1 

1 

1 

8.356031742 

10.70448409 

1 

-1 

-1 

1 

1 

1 

-1 

1 

0 

12.19686343 

1 

-1 

1 

-1 

-1 

-1 

-1 

1 

0.82080326 

10.67592725 

1 

-1 

1 

-1 

-1 

-1 

1 

-1 

0.220307227 

11.95604637 

-1 

1 

-1 

1 

1 

1 

-1 

1 

0 

12.27324157 

1 

-1 

-1 

-1 

1 

-1 

-1 

1 

0.82080326 

1.648543879 

-1 

-1 

1 

-1 

-1 

-1 

1 

1 

1.63713001 

10.19278707 

-1 

1 

1 

1 

-1 

1 

1 

-1 

0 

347.2821683 

1 

1 

1 

1 

-1 

1 

-1 

-1 

0 

348.9279031 

-1 

-1 

-1 

-1 

1 

-1 

1 

1 

12.85493611 

2.029011928 

-1 

-1 

-1 

1 

-1 

1 

1 

-1 

0 

13.89327317 

0 

0 

-1 

0 

0 

0 

0 

0 

0.004476521 

4.576830156 

-1 

-1 

-1 

1 

-1 

-1 

1 

1 

4.077927425 

6.010761982 

-1 

-1 

1 

1 

1 

1 

1 

-1 

0 

355.2176164 

1 

-1 

1 

1 

1 

1 

1 

1 

0 

283.1570151 

0 

0 

0 

0 

0 

0 

0 

1 

0.004476506 

8.650500413 

1 

-1 

-1 

-1 

-1 

-1 

-1 

-1 

0.820803265 

1.145741568 

0 

0 

0 

0 

0 

0 

0 

0 

0.004476511 

8.003210883 

1 

1 

-1 

-1 

-1 

-1 

-1 

1 

1.637130005 

1.113309233 

-1 

1 

1 

-1 

1 

-1 

-1 

-1 

1.63713001 

9.883301509 

-1 

1 

1 

-1 

-1 

1 

1 

1 

0 

12.58695421 

1 

1 

-1 

-1 

-1 

1 

1 

1 

0.004476511 

1.128043089 

0 

0 

0 

0 

0 

0 

-1 

0 

0.004476516 

18.16765092 

1 

1 

-1 

-1 

1 

1 

1 

-1 

0.004476511 

1.165215551 

1 

1 

-1 

-1 

1 

1 

-1 

1 

0.004476511 

1.147504641 

1 

1 

-1 

1 

-1 

1 

-1 

1 

0 

12.822396 

0 

-1 

0 

0 

0 

0 

0 

0 

0.004476506 

7.591606755 

-1 

-1 

-1 

-1 

-1 

-1 

1 

-1 

0.82080326 

1.111844522 

1 

-1 

1 

1 

-1 

-1 

1 

1 

0.312667506 

278.8550166 

-1 

-1 

1 

-1 

-1 

1 

1 

-1 

0 

13.66741485 

-1 

1 

1 

-1 

-1 

1 

-1 

-1 

0 

13.34852041 

-1 

1 

1 

-1 

-1 

-1 

-1 

1 

1.637130005 

8.078617813 

1 

-1 

1 

1 

-1 

1 

-1 

1 

0 

286.028231 

-1 

1 

-1 

-1 

-1 

-1 

-1 

-1 

0.820803255 

1.144141783 

1 

-1 

-1 

1 

-1 

1 

-1 

-1 

0 

13.61068865 

-1 

0 

0 

0 

0 

0 

0 

0 

20.50392892 

9.560084731 

1 

1 

1 

1 

1 

1 

-1 

1 

0 

267.8732819 

Table  7.2:  FCD  Experimental  Results  Part  2 
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alpha 

beta 

maxjterations 

num_ants 

num_neighbors 

q_0 

rho 

xi 

percent  over  optimal 

run  time 

-1 

1 

-1 

1 

-1 

1 

1 

1 

0 

11.55188327 

1 

-1 

1 

1 

1 

-1 

1 

-1 

0.004476511 

100.8008896 

-1 

-1 

-1 

1 

-1 

1 

-1 

1 

0 

11.57078539 

-1 

-1 

-1 

1 

1 

1 

-1 

-1 

0 

12.44630337 

1 

-1 

1 

-1 

1 

-1 

-1 

-1 

3.323575043 

10.61495162 

1 

1 

1 

1 

1 

-1 

-1 

-1 

8.356031737 

86.57042836 

-1 

-1 

1 

1 

1 

-1 

-1 

-1 

0.004476516 

95.10207429 

-1 

-1 

-1 

1 

1 

1 

1 

1 

0 

11.72501128 

-1 

-1 

1 

-1 

1 

1 

1 

1 

0 

11.98009328 

1 

1 

1 

-1 

1 

1 

1 

1 

0 

12.22406884 

-1 

1 

-1 

-1 

-1 

-1 

1 

1 

27.21905124 

1.097110806 

-1 

-1 

1 

-1 

1 

-1 

1 

-1 

36.7471159 

12.69746448 

-1 

-1 

1 

-1 

1 

1 

-1 

-1 

0 

12.84176397 

-1 

1 

1 

-1 

1 

1 

-1 

1 

0 

11.85530815 

1 

-1 

-1 

-1 

1 

-1 

1 

-1 

0.004476516 

1.869312161 

1 

1 

1 

-1 

-1 

1 

1 

-1 

0 

13.04782186 

1 

1 

1 

-1 

-1 

-1 

1 

1 

0.024175698 

13.86908252 

-1 

-1 

1 

1 

1 

1 

-1 

1 

0 

263.5471229 

-1 

-1 

-1 

1 

1 

-1 

1 

-1 

0.820803276 

14.00132335 

-1 

1 

1 

1 

1 

1 

1 

1 

0 

262.6201033 

-1 

1 

-1 

1 

1 

-1 

-1 

-1 

0.004476501 

11.84573556 

-1 

-1 

-1 

-1 

-1 

1 

1 

1 

0.004476511 

1.141610316 

-1 

1 

1 

1 

-1 

-1 

-1 

-1 

1.63713001 

274.9328212 

0 

0 

0 

1 

0 

0 

0 

0 

0.004476506 

39.06677134 

-1 

1 

-1 

1 

1 

1 

1 

-1 

0 

12.48501304 

1 

-1 

-1 

1 

-1 

-1 

-1 

1 

27.21905123 

6.271068669 

1 

-1 

-1 

-1 

-1 

1 

1 

-1 

0.004476511 

1.139471357 

-1 

1 

-1 

-1 

1 

1 

1 

1 

0.004476511 

1.10176461 

0 

0 

0 

0 

-1 

0 

0 

0 

4.959628171 

10.77465139 

1 

1 

1 

-1 

1 

-1 

1 

-1 

13.28078956 

9.88961111 

1 

1 

1 

-1 

1 

-1 

-1 

1 

14.50413045 

11.16779772 

1 

1 

-1 

1 

-1 

1 

1 

-1 

0 

12.52054144 

-1 

1 

-1 

-1 

-1 

1 

-1 

1 

0.004476511 

1.104824077 

1 

1 

-1 

1 

1 

1 

1 

1 

0 

11.57956528 

1 

-1 

1 

1 

-1 

-1 

-1 

-1 

4.077927425 

468.4784276 

1 

1 

-1 

-1 

1 

-1 

1 

1 

0.004476511 

1.96316543 

-1 

-1 

1 

1 

-1 

1 

-1 

-1 

0 

323.6814468 

-1 

-1 

-1 

-1 

-1 

1 

-1 

-1 

0.004476511 

1.102508 

1 

-1 

1 

-1 

-1 

1 

1 

1 

0 

11.88593299 

1 

1 

1 

1 

1 

1 

1 

-1 

0 

317.8588666 

-1 

-1 

-1 

-1 

1 

-1 

-1 

-1 

0.024175693 

1.711611252 

-1 

-1 

-1 

-1 

1 

1 

1 

-1 

0.004476511 

1.142997225 

0 

0 

0 

0 

0 

1 

0 

0 

0 

16.77154163 

0 

0 

0 

0 

0 

0 

0 

-1 

4.07792743 

8.818615355 

0 

1 

0 

0 

0 

0 

0 

0 

0.004476516 

14.23597804 

Table  7.3:  FCD  Experimental  Results  Part  3 
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